The aim of these exercises is to continue making sure you’re comfortable with handling multivariate data. In this chapter’s, you’ll focus on analyses based on dissimilarities/distances, including fitting linear models to these kinds of response variables.

For each of the examples below, you should follow the sequence we’ve used previously, as far as it’s sensible:

  1. What is the biological question?

  2. Is the predictor continuous or categorical?

  3. Write out the linear model corresponding to this question.

  4. What distribution do you expect the response variable to follow?

  5. What are the assumptions behind the statistical model you’ll fit?

    1. Are those assumptions satisfied?
  6. Fit the model

    1. How will you assess whether the model fits well?

    2. Can you detect an effect of the predictors?

    3. How do you measure the effect?

  7. What do you conclude (including any cautions)


A

Dixon et al. (2018) focused on assemblages of reptiles in forests and woodlands of southeastern Australia, with a particular interest in relationships between these assemblages and the fire history of the landscape. They identified 81 sites that varied in time since fire, from 6 months to >96 y. Rather than a continuum, three categories of time since fire were used, 0.5-2y, 6-12y, and >96y. Sites were also classified according to habitat.

Reptile assemblages were sampled with a range of methods, including visual surveys and camera traps, to give counts of 20 reptiles, whose abundance ranged from 0-1 to 0-126, depending on species.

Data are available from dryad, as Rep_abund.csv. You’ll want to focus on the columns with reptile numbers, and the one with the fire history (tsf). For convenience, we’ve extracted those data as dixonbiota.csv

Start by having a quick look at the data.

df <- read_csv("data/dixonbiota.csv")
Rows: 79 Columns: 22── Column specification ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr  (2): site, tsf
dbl (20): adup, aplat, amur, amac, aram, dcor, ecun, esax, eul, htal, ldel, lgui, lwhi, ppor, pent, pspe, ptex, rdie, tnig, vros
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(df, 10)

The file is pretty straightforward, with tsf being the fire category, and columns 3-22 each representing one reptile taxon.

Suppose you think of a reptile assemblage as simply the species that are present, regardless of how common they are. What would you conclude about the infuence of time since fire now?

Repeat the preceding analysis after creating a new data file that is presence-absence only.

df1.jac <- vegdist(df1,binary=TRUE,method='jaccard')
#Run 3d first
df1s.mds3 <- metaMDS(df1.jac,k=3,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.1404128 
Run 1 stress 0.1480196 
Run 2 stress 0.1405858 
... Procrustes: rmse 0.01565156  max resid 0.08607928 
Run 3 stress 0.1485349 
Run 4 stress 0.1417498 
Run 5 stress 0.1404104 
... New best solution
... Procrustes: rmse 0.0023304  max resid 0.01808758 
Run 6 stress 0.1404131 
... Procrustes: rmse 0.002334262  max resid 0.01834106 
Run 7 stress 0.1415289 
Run 8 stress 0.140411 
... Procrustes: rmse 0.0007835425  max resid 0.005200482 
... Similar to previous best
Run 9 stress 0.1486406 
Run 10 stress 0.1404103 
... New best solution
... Procrustes: rmse 0.001207819  max resid 0.009469643 
... Similar to previous best
Run 11 stress 0.1404117 
... Procrustes: rmse 0.000916645  max resid 0.006636852 
... Similar to previous best
Run 12 stress 0.1405581 
... Procrustes: rmse 0.01501831  max resid 0.08924602 
Run 13 stress 0.143846 
Run 14 stress 0.1495992 
Run 15 stress 0.140558 
... Procrustes: rmse 0.01458579  max resid 0.08787773 
Run 16 stress 0.1404148 
... Procrustes: rmse 0.0015373  max resid 0.01187026 
Run 17 stress 0.1448137 
Run 18 stress 0.1443774 
Run 19 stress 0.140411 
... Procrustes: rmse 0.001777414  max resid 0.01404256 
Run 20 stress 0.1500017 
Run 21 stress 0.1485371 
Run 22 stress 0.1454737 
Run 23 stress 0.140411 
... Procrustes: rmse 0.00189255  max resid 0.01470197 
Run 24 stress 0.1496201 
Run 25 stress 0.1405615 
... Procrustes: rmse 0.01445625  max resid 0.08690181 
Run 26 stress 0.1477962 
Run 27 stress 0.1478653 
Run 28 stress 0.1488089 
Run 29 stress 0.1444777 
Run 30 stress 0.1405647 
... Procrustes: rmse 0.01491986  max resid 0.08994996 
Run 31 stress 0.1405594 
... Procrustes: rmse 0.01486184  max resid 0.08799154 
Run 32 stress 0.1405612 
... Procrustes: rmse 0.01487556  max resid 0.0876917 
Run 33 stress 0.140559 
... Procrustes: rmse 0.01485206  max resid 0.0880716 
Run 34 stress 0.1425376 
Run 35 stress 0.140417 
... Procrustes: rmse 0.001866647  max resid 0.01406319 
Run 36 stress 0.1418405 
Run 37 stress 0.1422689 
Run 38 stress 0.1486242 
Run 39 stress 0.1413517 
Run 40 stress 0.1468941 
*** Best solution repeated 2 times
stressplot(df1s.mds3, main="Shepard plot")

df1s.mds3

Call:
metaMDS(comm = df1.jac, k = 3, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.jac 
Distance: binary jaccard 

Dimensions: 3 
Stress:     0.1404103 
Stress type 1, weak ties
Best solution was repeated 2 times in 40 tries
The best solution was from try 10 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing
# Now look at 2d
df1s.mds2 <- metaMDS(df1.jac,k=2,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.2009509 
Run 1 stress 0.2079375 
Run 2 stress 0.2224574 
Run 3 stress 0.2162073 
Run 4 stress 0.225869 
Run 5 stress 0.2284785 
Run 6 stress 0.2305984 
Run 7 stress 0.2199417 
Run 8 stress 0.202096 
Run 9 stress 0.2011626 
... Procrustes: rmse 0.006551419  max resid 0.03916879 
Run 10 stress 0.2104293 
Run 11 stress 0.2015109 
Run 12 stress 0.2009556 
... Procrustes: rmse 0.001703517  max resid 0.01335007 
Run 13 stress 0.2135392 
Run 14 stress 0.2246277 
Run 15 stress 0.2240908 
Run 16 stress 0.2250626 
Run 17 stress 0.200954 
... Procrustes: rmse 0.002998174  max resid 0.01914574 
Run 18 stress 0.2095741 
Run 19 stress 0.2186864 
Run 20 stress 0.203307 
Run 21 stress 0.2081138 
Run 22 stress 0.2212435 
Run 23 stress 0.2132584 
Run 24 stress 0.2055683 
Run 25 stress 0.2019164 
Run 26 stress 0.2108398 
Run 27 stress 0.2151842 
Run 28 stress 0.2009453 
... New best solution
... Procrustes: rmse 0.002748903  max resid 0.0197823 
Run 29 stress 0.206872 
Run 30 stress 0.2224539 
Run 31 stress 0.2031111 
Run 32 stress 0.2064081 
Run 33 stress 0.2149136 
Run 34 stress 0.2009451 
... New best solution
... Procrustes: rmse 0.002076934  max resid 0.0168433 
Run 35 stress 0.2111817 
Run 36 stress 0.2164562 
Run 37 stress 0.2027405 
Run 38 stress 0.2161414 
Run 39 stress 0.2242858 
Run 40 stress 0.2071548 
Run 41 stress 0.203393 
Run 42 stress 0.2266633 
Run 43 stress 0.2121352 
Run 44 stress 0.2175283 
Run 45 stress 0.2104325 
Run 46 stress 0.2029998 
Run 47 stress 0.2159185 
Run 48 stress 0.2042398 
Run 49 stress 0.2134569 
Run 50 stress 0.2145411 
Run 51 stress 0.2081612 
Run 52 stress 0.223412 
Run 53 stress 0.2327611 
Run 54 stress 0.2199861 
Run 55 stress 0.2115448 
Run 56 stress 0.209779 
Run 57 stress 0.2019152 
Run 58 stress 0.2236233 
Run 59 stress 0.2015106 
Run 60 stress 0.2036071 
Run 61 stress 0.2118656 
Run 62 stress 0.2120415 
Run 63 stress 0.2222671 
Run 64 stress 0.2176459 
Run 65 stress 0.2133805 
Run 66 stress 0.2009556 
... Procrustes: rmse 0.003462696  max resid 0.0201114 
Run 67 stress 0.201907 
Run 68 stress 0.2009456 
... Procrustes: rmse 0.0003121967  max resid 0.001856644 
... Similar to previous best
*** Best solution repeated 1 times
stressplot(df1s.mds2, main="Shepard plot")

df1s.mds2

Call:
metaMDS(comm = df1.jac, k = 2, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.jac 
Distance: binary jaccard 

Dimensions: 2 
Stress:     0.2009451 
Stress type 1, weak ties
Best solution was repeated 1 time in 68 tries
The best solution was from try 34 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing

Marginal decision here - stress just over 0.2 for k=2, so might be OK. Stress good for k=3, but 2-d MDS usually better for reporting or showing to audience

a<-as.data.frame(scores(df1s.mds3))
a<-cbind(df[c(1:2)],a)   #Add site names & symbols from original data file
p1<-ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=tsf, ) )+
  geom_point()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=sym3,
                     name="TSF",
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=a, aes(x=NMDS1, y=NMDS3, color = tsf, ) )+
  geom_point()+
  labs(y="MDS3", x="MDS1")+
  scale_shape_manual(values=sym3,
                     name="TSF",
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )

p1 + p2 + plot_layout(guides='collect') & theme_qk() & scale_color_uchicago()

Unburned sites separate very clearly on MDS1. Pattern clearer than with abundances included

Do permanova including tsf

df1.ado <- adonis2(df1.jac~tsf,data=df,permutations=999)
print(df1.ado)
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Permutation: free
Number of permutations: 999

adonis2(formula = df1.jac ~ tsf, data = df, permutations = 999)
         Df SumOfSqs      R2      F Pr(>F)    
tsf       2   3.7853 0.20678 9.9057  0.001 ***
Residual 76  14.5210 0.79322                  
Total    78  18.3062 1.00000                  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Clear fire signal

Extension activity

Dixon and her colleagues also recorded a range of habitat variables, which we’ve collected for you in dixonenv.csv

dixonenv <- read_csv("data/dixonenv.csv")
Rows: 79 Columns: 14── Column specification ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr  (4): site, tsf, veg, aspect
dbl (10): lcwd, shrcov, grcov, litcov, litdep, rocks, elev, warm, cold, twi
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(dixonenv, 10)

They recorded Coarse Woody Debris, which was long-transformed, litter cover (litcov), % cover of groundcover (grcov), and % cover of shrubs (shrcov), and rock cover. In the original paper, Dixon et al. were comfortable that these predictors weren’t correlated.

mvabund

dixonbiotamv <- mvabund(dixonbiota[,-(1:2)])
dixonbiotamv.mv <- manyglm(dixonbiotamv~dixonenv$tsf+dixonenv$lcwd+dixonenv$litcov+dixonenv$grcov+dixonenv$shrcov+dixonenv$lrocks,family="poisson")
plot(dixonbiotamv.mv)

# still hetero vars so use -ve binomial
dixonbiotamv1.mv <- manyglm(dixonbiotamv~dixonenv$tsf+dixonenv$lcwd+dixonenv$litcov+dixonenv$grcov+dixonenv$shrcov+dixonenv$lrocks,family="negative_binomial")
plot(dixonbiotamv1.mv)

anova(dixonbiotamv1.mv)
Time elapsed: 0 hr 1 min 39 sec
Analysis of Deviance Table

Model: dixonbiotamv ~ dixonenv$tsf + dixonenv$lcwd + dixonenv$litcov + dixonenv$grcov + dixonenv$shrcov + dixonenv$lrocks

Multivariate test:
                Res.Df Df.diff   Dev Pr(>Dev)    
(Intercept)         78                           
dixonenv$tsf        76       2 331.6    0.001 ***
dixonenv$lcwd       75       1  61.5    0.001 ***
dixonenv$litcov     74       1  25.9    0.265    
dixonenv$grcov      73       1  36.4    0.046 *  
dixonenv$shrcov     72       1  21.5    0.607    
dixonenv$lrocks     71       1  58.0    0.001 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Arguments:
 Test statistics calculated assuming uncorrelated response (for faster computation) 
 P-value calculated using 999 iterations via PIT-trap resampling.

B Simple MDS & Permanova

Let’s return to the Hutto and Barrett (2021) example from the Chapter 15 exercises. There, we used similarity-based analyses to explore the relationship between frog assemblages in ponds and the degree of urbanization surrounding. This seems a question that could just as easily be examined using distances or dissimilarities.

Return to the data and assess whether frog assemblages (using the standardized abundance scale or presence-absence) differ with urbanization.

Did your conclusions differ from those obtained using RDA?

df <- read_csv("data/huttoamph.csv")
Rows: 50 Columns: 14── Column specification ───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr  (1): type
dbl (13): Site, ACRE, BAME, BFOW, GCAR, HCIN, HVER, LCAT, LCLA, LPAL, LSPH, PCRU, PFER
ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
head(df,10)
#ponds 9 and 40 had no frogs. Remove for distance analysis
df <- df %>% 
  filter(Site!= 9 & Site !=40)
df

The explanation of the variables is in the previous chapter’s exercises.

df1 <- df[,-(1:2)]   #Stripped back file without labels, leaving only 
df1.bc <- vegdist(df1,'bray')
#Run 2d first
df1.mds2 <- metaMDS(df1.bc,k=2,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.196708 
Run 1 stress 0.1877626 
... New best solution
... Procrustes: rmse 0.06205021  max resid 0.2289634 
Run 2 stress 0.1888251 
Run 3 stress 0.1920872 
Run 4 stress 0.1853642 
... New best solution
... Procrustes: rmse 0.07936891  max resid 0.4999309 
Run 5 stress 0.1855846 
... Procrustes: rmse 0.01953946  max resid 0.1056279 
Run 6 stress 0.1984394 
Run 7 stress 0.1862238 
Run 8 stress 0.200767 
Run 9 stress 0.1961235 
Run 10 stress 0.1853641 
... New best solution
... Procrustes: rmse 9.823368e-05  max resid 0.0004897303 
... Similar to previous best
Run 11 stress 0.204919 
Run 12 stress 0.2001183 
Run 13 stress 0.1875595 
Run 14 stress 0.2031338 
Run 15 stress 0.1875596 
Run 16 stress 0.2001202 
Run 17 stress 0.1877626 
Run 18 stress 0.1968146 
Run 19 stress 0.2191092 
Run 20 stress 0.1862237 
Run 21 stress 0.1862237 
Run 22 stress 0.2016092 
Run 23 stress 0.1919379 
Run 24 stress 0.1921465 
Run 25 stress 0.1875596 
Run 26 stress 0.1862241 
Run 27 stress 0.1888273 
Run 28 stress 0.1958875 
Run 29 stress 0.2103827 
Run 30 stress 0.196403 
Run 31 stress 0.1920872 
Run 32 stress 0.1853644 
... Procrustes: rmse 0.0001814753  max resid 0.0007769668 
... Similar to previous best
Run 33 stress 0.1929083 
Run 34 stress 0.1948815 
Run 35 stress 0.2023845 
Run 36 stress 0.1961143 
Run 37 stress 0.1978435 
Run 38 stress 0.214542 
Run 39 stress 0.1970125 
Run 40 stress 0.2055679 
*** Best solution repeated 2 times
stressplot(df1.mds2, main="Shepard plot")

df1.mds2

Call:
metaMDS(comm = df1.bc, k = 2, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.bc 
Distance: bray 

Dimensions: 2 
Stress:     0.1853641 
Stress type 1, weak ties
Best solution was repeated 2 times in 40 tries
The best solution was from try 10 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing
# Now look at 3d
df.mds3 <- metaMDS(df1.bc,k=3,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.129187 
Run 1 stress 0.1388768 
Run 2 stress 0.1382418 
Run 3 stress 0.1383885 
Run 4 stress 0.129186 
... New best solution
... Procrustes: rmse 0.001355541  max resid 0.007232405 
... Similar to previous best
Run 5 stress 0.129186 
... Procrustes: rmse 0.0005937864  max resid 0.001870774 
... Similar to previous best
Run 6 stress 0.1291858 
... New best solution
... Procrustes: rmse 0.0008376747  max resid 0.004337427 
... Similar to previous best
Run 7 stress 0.129186 
... Procrustes: rmse 0.0007498725  max resid 0.00420023 
... Similar to previous best
Run 8 stress 0.1291864 
... Procrustes: rmse 0.0002763668  max resid 0.0009518002 
... Similar to previous best
Run 9 stress 0.1383896 
Run 10 stress 0.1384317 
Run 11 stress 0.1289614 
... New best solution
... Procrustes: rmse 0.0109308  max resid 0.05605748 
Run 12 stress 0.128962 
... Procrustes: rmse 0.0009039712  max resid 0.005251417 
... Similar to previous best
Run 13 stress 0.1383929 
Run 14 stress 0.1402552 
Run 15 stress 0.128962 
... Procrustes: rmse 0.0009246251  max resid 0.005329771 
... Similar to previous best
Run 16 stress 0.1291859 
... Procrustes: rmse 0.01058901  max resid 0.0547769 
Run 17 stress 0.1289615 
... Procrustes: rmse 0.000419003  max resid 0.002034791 
... Similar to previous best
Run 18 stress 0.1380889 
Run 19 stress 0.129186 
... Procrustes: rmse 0.01092673  max resid 0.05630155 
Run 20 stress 0.1383886 
Run 21 stress 0.1289617 
... Procrustes: rmse 0.0002465603  max resid 0.0007331739 
... Similar to previous best
Run 22 stress 0.1291604 
... Procrustes: rmse 0.009890855  max resid 0.05184425 
Run 23 stress 0.1291601 
... Procrustes: rmse 0.009490803  max resid 0.05078321 
Run 24 stress 0.1291855 
... Procrustes: rmse 0.01066625  max resid 0.05496251 
Run 25 stress 0.1391225 
Run 26 stress 0.129186 
... Procrustes: rmse 0.01099671  max resid 0.05629543 
Run 27 stress 0.1291866 
... Procrustes: rmse 0.01114331  max resid 0.05711187 
Run 28 stress 0.1384282 
Run 29 stress 0.1291861 
... Procrustes: rmse 0.01101164  max resid 0.0563452 
Run 30 stress 0.1383979 
Run 31 stress 0.1289613 
... New best solution
... Procrustes: rmse 0.0004223583  max resid 0.001654659 
... Similar to previous best
Run 32 stress 0.1380874 
Run 33 stress 0.1291858 
... Procrustes: rmse 0.01077276  max resid 0.05558138 
Run 34 stress 0.1291859 
... Procrustes: rmse 0.01081878  max resid 0.05582879 
Run 35 stress 0.1291858 
... Procrustes: rmse 0.01042548  max resid 0.05433355 
Run 36 stress 0.1385874 
Run 37 stress 0.1291863 
... Procrustes: rmse 0.01092361  max resid 0.05638233 
Run 38 stress 0.1308828 
Run 39 stress 0.1384275 
Run 40 stress 0.1388191 
*** Best solution repeated 1 times
stressplot(df.mds3, main="Shepard plot")

df.mds3

Call:
metaMDS(comm = df1.bc, k = 3, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.bc 
Distance: bray 

Dimensions: 3 
Stress:     0.1289613 
Stress type 1, weak ties
Best solution was repeated 1 time in 40 tries
The best solution was from try 31 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing

2-d is acceptable

a<-as.data.frame(scores(df1.mds2))
a<-cbind(df[c(1:2)],a)   #Add site names & symbols from original data file
ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=type, ) )+
  geom_point()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=sym3,
                     name="Urbanization",
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )+
  theme_qk() + scale_color_uchicago()

df1.ado <- adonis2(df1.bc~type,data=df,permutations=999)
print(df1.ado)
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Permutation: free
Number of permutations: 999

adonis2(formula = df1.bc ~ type, data = df, permutations = 999)
         Df SumOfSqs      R2      F Pr(>F)  
type      2    0.912 0.07281 1.7669  0.045 *
Residual 45   11.614 0.92719                
Total    47   12.526 1.00000                
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Permanova detects a difference; MDS plot suggests low urbanization may be different from the other two, though worth looking at dispersion as well, as spread of this type is less than the other two on plot

Repeat analysis with presence-absence

df1.jac <- vegdist(df1,binary=TRUE,method=‘jaccard’)

#Run 2d first
df1.jac <- vegdist(df1,binary=TRUE,method='jaccard')
df1.mds2 <- metaMDS(df1.jac,k=2,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.1551612 
Run 1 stress 0.1630921 
Run 2 stress 0.1622263 
Run 3 stress 0.1724054 
Run 4 stress 0.1585799 
Run 5 stress 0.1761433 
Run 6 stress 0.1608258 
Run 7 stress 0.1734443 
Run 8 stress 0.1621126 
Run 9 stress 0.1536788 
... New best solution
... Procrustes: rmse 0.01901219  max resid 0.06506965 
Run 10 stress 0.1760283 
Run 11 stress 0.1609339 
Run 12 stress 0.1565712 
Run 13 stress 0.1939366 
Run 14 stress 0.1843833 
Run 15 stress 0.162526 
Run 16 stress 0.1622263 
Run 17 stress 0.1570453 
Run 18 stress 0.1708675 
Run 19 stress 0.1537515 
... Procrustes: rmse 0.003076067  max resid 0.01474863 
Run 20 stress 0.1539129 
... Procrustes: rmse 0.00947995  max resid 0.04425278 
Run 21 stress 0.1664178 
Run 22 stress 0.1795028 
Run 23 stress 0.1585801 
Run 24 stress 0.1555104 
Run 25 stress 0.1708671 
Run 26 stress 0.167874 
Run 27 stress 0.1673092 
Run 28 stress 0.1719174 
Run 29 stress 0.1591019 
Run 30 stress 0.1711482 
Run 31 stress 0.1795847 
Run 32 stress 0.1608763 
Run 33 stress 0.1721431 
Run 34 stress 0.1585798 
Run 35 stress 0.1850134 
Run 36 stress 0.1607836 
Run 37 stress 0.1586731 
Run 38 stress 0.1738412 
Run 39 stress 0.2004114 
Run 40 stress 0.1672269 
Run 41 stress 0.1585799 
Run 42 stress 0.1591019 
Run 43 stress 0.162564 
Run 44 stress 0.1620556 
Run 45 stress 0.1679859 
Run 46 stress 0.171033 
Run 47 stress 0.1617404 
Run 48 stress 0.1538246 
... Procrustes: rmse 0.01103138  max resid 0.05221548 
Run 49 stress 0.172142 
Run 50 stress 0.1960733 
Run 51 stress 0.1731862 
Run 52 stress 0.1627813 
Run 53 stress 0.1536788 
... Procrustes: rmse 0.0001144048  max resid 0.000322326 
... Similar to previous best
*** Best solution repeated 1 times
stressplot(df1.mds2, main="Shepard plot")

df1.mds2

Call:
metaMDS(comm = df1.jac, k = 2, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.jac 
Distance: binary jaccard 

Dimensions: 2 
Stress:     0.1536788 
Stress type 1, weak ties
Best solution was repeated 1 time in 53 tries
The best solution was from try 9 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing
# Now look at 3d
df.mds3 <- metaMDS(df1.jac,k=3,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.09466813 
Run 1 stress 0.09466802 
... New best solution
... Procrustes: rmse 0.0004764601  max resid 0.002102695 
... Similar to previous best
Run 2 stress 0.09466791 
... New best solution
... Procrustes: rmse 0.0001518385  max resid 0.0004555261 
... Similar to previous best
Run 3 stress 0.09466803 
... Procrustes: rmse 0.0001108676  max resid 0.0004745956 
... Similar to previous best
Run 4 stress 0.1005052 
Run 5 stress 0.09466787 
... New best solution
... Procrustes: rmse 5.340853e-05  max resid 0.0001992126 
... Similar to previous best
Run 6 stress 0.1005044 
Run 7 stress 0.1005056 
Run 8 stress 0.09466803 
... Procrustes: rmse 0.0001288817  max resid 0.0005023949 
... Similar to previous best
Run 9 stress 0.1005047 
Run 10 stress 0.09466793 
... Procrustes: rmse 5.310415e-05  max resid 0.0002709636 
... Similar to previous best
Run 11 stress 0.09466796 
... Procrustes: rmse 0.0001256892  max resid 0.0003435652 
... Similar to previous best
Run 12 stress 0.094668 
... Procrustes: rmse 9.745181e-05  max resid 0.0004669281 
... Similar to previous best
Run 13 stress 0.09693819 
Run 14 stress 0.09694632 
Run 15 stress 0.1005051 
Run 16 stress 0.09466833 
... Procrustes: rmse 0.0002133226  max resid 0.001028399 
... Similar to previous best
Run 17 stress 0.09466821 
... Procrustes: rmse 0.0002054658  max resid 0.0009036094 
... Similar to previous best
Run 18 stress 0.09466811 
... Procrustes: rmse 0.0001423043  max resid 0.0004843394 
... Similar to previous best
Run 19 stress 0.09466805 
... Procrustes: rmse 0.0003941737  max resid 0.002017784 
... Similar to previous best
Run 20 stress 0.109714 
Run 21 stress 0.1005051 
Run 22 stress 0.09466816 
... Procrustes: rmse 0.000195071  max resid 0.0006552548 
... Similar to previous best
Run 23 stress 0.09466783 
... New best solution
... Procrustes: rmse 6.655458e-05  max resid 0.0002929959 
... Similar to previous best
Run 24 stress 0.09466805 
... Procrustes: rmse 0.0001689595  max resid 0.0006247545 
... Similar to previous best
Run 25 stress 0.1005048 
Run 26 stress 0.09466801 
... Procrustes: rmse 0.0001773767  max resid 0.0009243069 
... Similar to previous best
Run 27 stress 0.09694616 
Run 28 stress 0.1004898 
Run 29 stress 0.09466783 
... New best solution
... Procrustes: rmse 0.0001519669  max resid 0.0004171994 
... Similar to previous best
Run 30 stress 0.1005053 
Run 31 stress 0.1097146 
Run 32 stress 0.1004899 
Run 33 stress 0.09466786 
... Procrustes: rmse 0.0001828158  max resid 0.000484661 
... Similar to previous best
Run 34 stress 0.1090648 
Run 35 stress 0.09466812 
... Procrustes: rmse 0.0002666448  max resid 0.0009620816 
... Similar to previous best
Run 36 stress 0.09694549 
Run 37 stress 0.09694624 
Run 38 stress 0.094668 
... Procrustes: rmse 0.0002647032  max resid 0.001383475 
... Similar to previous best
Run 39 stress 0.09466798 
... Procrustes: rmse 0.0002442504  max resid 0.0007524015 
... Similar to previous best
Run 40 stress 0.09466802 
... Procrustes: rmse 0.000254564  max resid 0.0008115543 
... Similar to previous best
*** Best solution repeated 6 times
stressplot(df.mds3, main="Shepard plot")

df.mds3

Call:
metaMDS(comm = df1.jac, k = 3, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1.jac 
Distance: binary jaccard 

Dimensions: 3 
Stress:     0.09466783 
Stress type 1, weak ties
Best solution was repeated 6 times in 40 tries
The best solution was from try 29 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing

Again, stress OK with 2D, good with 3D. Look at 2D first

a<-as.data.frame(scores(df1.mds2))
a<-cbind(df[c(1:2)],a)   #Add site names & symbols from original data file
ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=type, ) )+
  geom_point()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=sym3,
                     name="Urbanization",
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )+
  theme_qk() + scale_color_uchicago()

df1.ado <- adonis2(df1.jac~type,data=df,permutations=999)
print(df1.ado)
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Permutation: free
Number of permutations: 999

adonis2(formula = df1.jac ~ type, data = df, permutations = 999)
         Df SumOfSqs      R2      F Pr(>F)
type      2   0.4525 0.04389 1.0328  0.403
Residual 45   9.8587 0.95611              
Total    47  10.3112 1.00000              

No separation seen on Presence-Absence

C

Griffen et al. (2019) examined the oral microbiome of humans, with a focus on the effects of HIV and AntiRetroviral therapy (ART) on this microbiome. They described the microbiome of 341 patients, who fell into one of two categories: HIV-, and HIV+ with ART. They also matched their samples as far as possible for sex, along with other conditions that can influence oral microbiomes (Candida infection, current smoking). We’ll use their records for these other categories as well. There were more HIV+ than - subjects, but within these groups, approximately equal numbers of two sexes, two Candida infection status, and two smoking categories.

The bacteriome was recorded separately using 16S RNA sequencing, which overed just over 600 taxa.

The data are available from their paper, as two Excel files in the supplementary information. The metadata gives HIV status, and a raft of demographic information. We’ll just work with HIV, sex, Candida, and current smoking. A second sheet has the bacteriome. The code chunk below reads this file (from a Downloads folder - you’ll need to modify the file location). It also screens the file for any bacterial taxa that were not present in this particular study, which reduces the bacterial taxa to 599.

library(readxl)
#Read in metadata. Lots of information we don't need, so a couple of iterations to get down to a few columns we want - HIV status, candida, current smoking, gender
df <- read_excel("~/Downloads/41598_2019_55703_MOESM3_ESM.xlsx", 
     range = "A1:S342", na = "NA")
df<-select(df, HIV, candida, gender, smoking_current)
head(df)
#df1 is the bacterial data
df1 <- read_excel("~/Downloads/41598_2019_55703_MOESM2_ESM.xlsx")
New names:
head(df1,10)
df1 <- df1 %>%
  select(where( ~ is.numeric(.x) && sum(.x) != 0))  #Drops any bacterial taxon not recorded in this study (i.e. where it's all zeroes)
df1 <- df1[,-1] #remove col1, which is sample ID

Use a dissimilarity based approach to assess the combined effects of HIV-ART, sex, Candida infection, and current smoking on the bacteriome

Do these factors act independently on the bacteriome?

Which effects are largest (use some graphical methods)?

First steps: Think about standardization and which distance measure you’ll use. Bray-Curtis seems common for these kind of data, and and counts of different taxa vary widely.

df1s <- wisconsin(df1)   #Common to do some standardisation, and counts of different taxa vary widely
df1s.bc <- vegdist(df1s,'bray')   #With standardisation (wisconsin)

Run factorial model using permanova, based on B-C

df1.ado <- adonis2(df1s.bc~HIV*candida*gender*smoking_current,data=df,permutations=999)
print(df1.ado)
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Permutation: free
Number of permutations: 999

adonis2(formula = df1s.bc ~ HIV * candida * gender * smoking_current, data = df, permutations = 999)
                                    Df SumOfSqs      R2      F Pr(>F)    
HIV                                  1    1.369 0.01250 4.3776  0.001 ***
candida                              1    1.461 0.01333 4.6694  0.001 ***
gender                               1    0.492 0.00449 1.5737  0.009 ** 
smoking_current                      1    1.111 0.01014 3.5513  0.001 ***
HIV:candida                          1    0.318 0.00290 1.0151  0.401    
HIV:gender                           1    0.354 0.00324 1.1332  0.204    
candida:gender                       1    0.257 0.00234 0.8202  0.861    
HIV:smoking_current                  1    0.269 0.00245 0.8584  0.765    
candida:smoking_current              1    0.344 0.00314 1.0989  0.243    
gender:smoking_current               1    0.394 0.00360 1.2602  0.098 .  
HIV:candida:gender                   1    0.349 0.00319 1.1164  0.203    
HIV:candida:smoking_current          1    0.280 0.00255 0.8943  0.707    
HIV:gender:smoking_current           1    0.265 0.00242 0.8463  0.793    
candida:gender:smoking_current       1    0.254 0.00232 0.8110  0.841    
HIV:candida:gender:smoking_current   1    0.362 0.00330 1.1562  0.213    
Residual                           325  101.664 0.92809                  
Total                              340  109.541 1.00000                  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Model shows main effects, but no interactions important (or at least, detected)

Could you fit a simpler model to the data?

Run simpler model

df2.ado <- adonis2(df1s.bc~HIV+candida+gender+smoking_current,data=df,permutations=999)
print(df2.ado)
Permutation test for adonis under reduced model
Terms added sequentially (first to last)
Permutation: free
Number of permutations: 999

adonis2(formula = df1s.bc ~ HIV + candida + gender + smoking_current, data = df, permutations = 999)
                 Df SumOfSqs      R2      F Pr(>F)    
HIV               1    1.369 0.01250 4.3775  0.001 ***
candida           1    1.461 0.01333 4.6693  0.001 ***
gender            1    0.492 0.00449 1.5736  0.013 *  
smoking_current   1    1.111 0.01014 3.5512  0.001 ***
Residual        336  105.108 0.95953                  
Total           340  109.541 1.00000                  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

No change; not too surprising, as all four main effects were detected easily.

Now go to data visualization

#Run 3d first
df1s.mds3 <- metaMDS(df1s.bc,k=3,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.1573956 
Run 1 stress 0.1579391 
Run 2 stress 0.1572758 
... New best solution
... Procrustes: rmse 0.01025549  max resid 0.1525195 
Run 3 stress 0.1579179 
Run 4 stress 0.157493 
... Procrustes: rmse 0.007283638  max resid 0.1009598 
Run 5 stress 0.1578733 
Run 6 stress 0.1593345 
Run 7 stress 0.1596534 
Run 8 stress 0.1578064 
Run 9 stress 0.1576408 
... Procrustes: rmse 0.005563033  max resid 0.08072257 
Run 10 stress 0.1574589 
... Procrustes: rmse 0.007053274  max resid 0.1007079 
Run 11 stress 0.1574113 
... Procrustes: rmse 0.008390411  max resid 0.1530435 
Run 12 stress 0.1584759 
Run 13 stress 0.1590917 
Run 14 stress 0.1579313 
Run 15 stress 0.1577686 
... Procrustes: rmse 0.00731071  max resid 0.08902012 
Run 16 stress 0.1579681 
Run 17 stress 0.157889 
Run 18 stress 0.1576122 
... Procrustes: rmse 0.01046473  max resid 0.1079021 
Run 19 stress 0.1575193 
... Procrustes: rmse 0.01344117  max resid 0.1517403 
Run 20 stress 0.1590183 
Run 21 stress 0.1572873 
... Procrustes: rmse 0.005938186  max resid 0.1084724 
Run 22 stress 0.1579056 
Run 23 stress 0.1578609 
Run 24 stress 0.1575202 
... Procrustes: rmse 0.01348549  max resid 0.1517058 
Run 25 stress 0.1580688 
Run 26 stress 0.1573537 
... Procrustes: rmse 0.005408552  max resid 0.09851134 
Run 27 stress 0.1584273 
Run 28 stress 0.1575528 
... Procrustes: rmse 0.007336861  max resid 0.1083012 
Run 29 stress 0.1591385 
Run 30 stress 0.1575529 
... Procrustes: rmse 0.007325558  max resid 0.1082856 
Run 31 stress 0.1584616 
Run 32 stress 0.157776 
Run 33 stress 0.1579148 
Run 34 stress 0.1590267 
Run 35 stress 0.1574802 
... Procrustes: rmse 0.01063831  max resid 0.1527981 
Run 36 stress 0.1574825 
... Procrustes: rmse 0.0125045  max resid 0.1520643 
Run 37 stress 0.158388 
Run 38 stress 0.1581063 
Run 39 stress 0.1586636 
Run 40 stress 0.159344 
Run 41 stress 0.1575585 
... Procrustes: rmse 0.007481761  max resid 0.1081395 
Run 42 stress 0.1580441 
Run 43 stress 0.1583809 
Run 44 stress 0.1586431 
Run 45 stress 0.1575179 
... Procrustes: rmse 0.009346751  max resid 0.1082067 
Run 46 stress 0.1585796 
Run 47 stress 0.1578686 
Run 48 stress 0.1595925 
Run 49 stress 0.1579636 
Run 50 stress 0.1572746 
... New best solution
... Procrustes: rmse 0.0005494158  max resid 0.00824335 
... Similar to previous best
*** Best solution repeated 1 times
stressplot(df1s.mds3, main="Shepard plot")

df1s.mds3

Call:
metaMDS(comm = df1s.bc, k = 3, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1s.bc 
Distance: bray 

Dimensions: 3 
Stress:     0.1572746 
Stress type 1, weak ties
Best solution was repeated 1 time in 50 tries
The best solution was from try 50 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing
# Now look at 2d
df1s.mds2 <- metaMDS(df1s.bc,k=2,autotransform=FALSE,try=40,trymax=80,maxit=200)
Run 0 stress 0.2249552 
Run 1 stress 0.2254408 
... Procrustes: rmse 0.01613467  max resid 0.2308209 
Run 2 stress 0.2253412 
... Procrustes: rmse 0.01606138  max resid 0.2302704 
Run 3 stress 0.2250828 
... Procrustes: rmse 0.009539203  max resid 0.1512118 
Run 4 stress 0.2250839 
... Procrustes: rmse 0.009650024  max resid 0.1525654 
Run 5 stress 0.2250832 
... Procrustes: rmse 0.009612515  max resid 0.1518641 
Run 6 stress 0.2248224 
... New best solution
... Procrustes: rmse 0.00631171  max resid 0.08263297 
Run 7 stress 0.2253967 
Run 8 stress 0.225527 
Run 9 stress 0.2248128 
... New best solution
... Procrustes: rmse 0.0004834727  max resid 0.004762261 
... Similar to previous best
Run 10 stress 0.2255447 
Run 11 stress 0.2249136 
... Procrustes: rmse 0.00484646  max resid 0.08193293 
Run 12 stress 0.2500137 
Run 13 stress 0.2256121 
Run 14 stress 0.2254928 
Run 15 stress 0.2258491 
Run 16 stress 0.2248191 
... Procrustes: rmse 0.0008375058  max resid 0.008573354 
... Similar to previous best
Run 17 stress 0.225849 
Run 18 stress 0.225634 
Run 19 stress 0.2249543 
... Procrustes: rmse 0.006228313  max resid 0.08122835 
Run 20 stress 0.2250713 
... Procrustes: rmse 0.008311559  max resid 0.1485133 
Run 21 stress 0.2258052 
Run 22 stress 0.2248935 
... Procrustes: rmse 0.004390114  max resid 0.08026974 
Run 23 stress 0.2257226 
Run 24 stress 0.2253296 
Run 25 stress 0.2257815 
Run 26 stress 0.2254095 
Run 27 stress 0.2257822 
Run 28 stress 0.2248496 
... Procrustes: rmse 0.004509314  max resid 0.0818388 
Run 29 stress 0.2250698 
... Procrustes: rmse 0.008236929  max resid 0.1489688 
Run 30 stress 0.225152 
... Procrustes: rmse 0.01050374  max resid 0.1503719 
Run 31 stress 0.2426061 
Run 32 stress 0.2253598 
Run 33 stress 0.418577 
Run 34 stress 0.2254408 
Run 35 stress 0.2254011 
Run 36 stress 0.2257955 
Run 37 stress 0.2256214 
Run 38 stress 0.22512 
... Procrustes: rmse 0.009424238  max resid 0.1483398 
Run 39 stress 0.2255057 
Run 40 stress 0.225153 
... Procrustes: rmse 0.0106079  max resid 0.1526557 
*** Best solution repeated 2 times
stressplot(df1s.mds2, main="Shepard plot")

df1s.mds2

Call:
metaMDS(comm = df1s.bc, k = 2, try = 40, trymax = 80, autotransform = FALSE,      maxit = 200) 

global Multidimensional Scaling using monoMDS

Data:     df1s.bc 
Distance: bray 

Dimensions: 2 
Stress:     0.2248128 
Stress type 1, weak ties
Best solution was repeated 2 times in 40 tries
The best solution was from try 9 (random start)
Scaling: centring, PC rotation, halfchange scaling 
Species: scores missing

3D model fits better than 2D. Stress for 2D above 0.2

Start with plots where symbol colour indicates levels of one of the factors

a<-as.data.frame(scores(df1s.mds3))
a<-cbind(df[c(1:4)],a)   #Add site names & symbols from original data file
p1<-ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=HIV) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=a, aes(x=NMDS1, y=NMDS3, color = HIV) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS3", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )

p1 + p2 + plot_layout(guides='collect') + plot_annotation(title="HIV status") & theme_qk() & scale_color_uchicago()


p1<-ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=candida) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=a, aes(x=NMDS1, y=NMDS3, color = candida) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS3", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )

p1 + p2 + plot_layout(guides='collect') + plot_annotation(title="Candida") & theme_qk() & scale_color_uchicago()


p1<-ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=gender) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=a, aes(x=NMDS1, y=NMDS3, color = gender) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS3", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )

p1 + p2 + plot_layout(guides='collect') + plot_annotation(title="Sex") & theme_qk() & scale_color_uchicago()


p1<-ggplot(data=a, aes(x=NMDS1, y=NMDS2, color=smoking_current) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=a, aes(x=NMDS1, y=NMDS3, color = smoking_current) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS3", x="MDS1")+
  scale_shape_manual(values=c(16,17),
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )

p1 + p2 + plot_layout(guides='collect') + plot_annotation(title="Current smoking") & theme_qk() & scale_color_uchicago()

Try some visualizations of pairs of factors using MDS1 & MDS2 (while 3D is preferred, stress of 3D is around .15, vs 0.22, so not huge difference) Separate HIV+/-, then use symbol colour to separate second factor.

p1<-ggplot(data=subset(a, candida=="Yes"), aes(x=NMDS1, y=NMDS2, color=HIV) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title="Candida")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=subset(a, candida =="No"), aes(x=NMDS1, y=NMDS2, color = HIV) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title = "No Candida")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )
p_combined <- p1 + p2 + plot_layout(guides='collect') & theme_qk() & scale_color_uchicago()
p_ranges_x <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_x[[1]]$range$range,
  ggplot_build(p_combined[[2]])$layout$panel_scales_x[[1]]$range$range)

p_ranges_y <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_y[[1]]$range$range,
                ggplot_build(p_combined[[2]])$layout$panel_scales_y[[1]]$range$range)
p_combined & 
  xlim(min(p_ranges_x), max(p_ranges_x)) & 
  ylim(min(p_ranges_y), max(p_ranges_y))

p1<-ggplot(data=subset(a, HIV=="Yes"), aes(x=NMDS1, y=NMDS2, color=candida) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title="HIV+")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=subset(a, HIV =="No"), aes(x=NMDS1, y=NMDS2, color = candida) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title = "HIV-")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )
p_combined <- p1 + p2 + plot_layout(guides='collect') & theme_qk() & scale_color_uchicago()
p_ranges_x <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_x[[1]]$range$range,
  ggplot_build(p_combined[[2]])$layout$panel_scales_x[[1]]$range$range)

p_ranges_y <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_y[[1]]$range$range,
                ggplot_build(p_combined[[2]])$layout$panel_scales_y[[1]]$range$range)
p_combined & 
  xlim(min(p_ranges_x), max(p_ranges_x)) & 
  ylim(min(p_ranges_y), max(p_ranges_y))

p1<-ggplot(data=subset(a, HIV=="Yes"), aes(x=NMDS1, y=NMDS2, color = smoking_current) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title="HIV+")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6), 
                                      title=NULL)
                     )
p2<-ggplot(data=subset(a, HIV =="No"), aes(x=NMDS1, y=NMDS2, color = smoking_current) )+
  geom_point()+
  stat_ellipse()+
  labs(y="MDS2", x="MDS1",
       title = "HIV-")+
  scale_shape_manual(
                     guide =
                         guide_legend(label.theme = element_text(size=6),
                                    title=NULL)
                     )
p_combined <- p1 + p2 + plot_layout(guides='collect') & theme_qk() & scale_color_uchicago()
p_ranges_x <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_x[[1]]$range$range,
  ggplot_build(p_combined[[2]])$layout$panel_scales_x[[1]]$range$range)

p_ranges_y <- c(ggplot_build(p_combined[[1]])$layout$panel_scales_y[[1]]$range$range,
                ggplot_build(p_combined[[2]])$layout$panel_scales_y[[1]]$range$range)
p_combined & 
  xlim(min(p_ranges_x), max(p_ranges_x)) & 
  ylim(min(p_ranges_y), max(p_ranges_y))

References

Dixon, Kelly M., Geoffrey J. Cary, Graeme L. Worboys, and Philip Gibbons. 2018. “The Disproportionate Importance of Long‐unburned Forests and Woodlands for Reptiles.” Ecology and Evolution 8 (22): 10952–63. https://doi.org/gfs587.
Griffen, Ann L., Zachary A. Thompson, Clifford J. Beall, Elizabeth A. Lilly, Carolina Granada, Kelly D. Treas, Kenneth R. DuBois, et al. 2019. “Significant Effect of HIV/HAART on Oral Microbiota Using Multivariate Analysis.” Scientific Reports 9 (1): 19946. https://doi.org/gsfpz3.
Hutto, David, and Kyle Barrett. 2021. “Do Urban Open Spaces Provide Refugia for Frogs in Urban Environments?” Plos One 16 (1). https://doi.org/gr2r9n.
LS0tCnRpdGxlOiAiQ2ggMTYgZXhlcmNpc2VzIgpvdXRwdXQ6IAogIGh0bWxfbm90ZWJvb2s6CiAgICB0aGVtZTogZmxhdGx5CmJpYmxpb2dyYXBoeTogd2ViX2V4LmJpYgotLS0KCmBgYHtyIGVjaG89RkFMU0UsIHJlc3VsdHM9J2hpZGUnfQojIFBhY2thZ2VzOiB2ZWdhbixtdmFidW5kCnNvdXJjZSgiUi9saWJyYXJpZXMuUiIpCnNvdXJjZSgiUi9hcHBlYXJhbmNlLlIiKQpsaWJyYXJ5KHZlZ2FuKQpsaWJyYXJ5KG12YWJ1bmQpCmBgYAoKVGhlIGFpbSBvZiB0aGVzZSBleGVyY2lzZXMgaXMgdG8gY29udGludWUgbWFraW5nIHN1cmUgeW91J3JlIGNvbWZvcnRhYmxlIHdpdGggaGFuZGxpbmcgbXVsdGl2YXJpYXRlIGRhdGEuIEluIHRoaXMgY2hhcHRlcidzLCB5b3UnbGwgZm9jdXMgb24gYW5hbHlzZXMgYmFzZWQgb24gZGlzc2ltaWxhcml0aWVzL2Rpc3RhbmNlcywgaW5jbHVkaW5nIGZpdHRpbmcgbGluZWFyIG1vZGVscyB0byB0aGVzZSBraW5kcyBvZiByZXNwb25zZSB2YXJpYWJsZXMuCgpGb3IgZWFjaCBvZiB0aGUgZXhhbXBsZXMgYmVsb3csIHlvdSBzaG91bGQgZm9sbG93IHRoZSBzZXF1ZW5jZSB3ZSd2ZSB1c2VkIHByZXZpb3VzbHksIGFzIGZhciBhcyBpdCdzIHNlbnNpYmxlOgoKMS4gIFdoYXQgaXMgdGhlIGJpb2xvZ2ljYWwgcXVlc3Rpb24/CgoyLiAgSXMgdGhlIHByZWRpY3RvciBjb250aW51b3VzIG9yIGNhdGVnb3JpY2FsPwoKMy4gIFdyaXRlIG91dCB0aGUgbGluZWFyIG1vZGVsIGNvcnJlc3BvbmRpbmcgdG8gdGhpcyBxdWVzdGlvbi4KCjQuICBXaGF0IGRpc3RyaWJ1dGlvbiBkbyB5b3UgZXhwZWN0IHRoZSByZXNwb25zZSB2YXJpYWJsZSB0byBmb2xsb3c/Cgo1LiAgV2hhdCBhcmUgdGhlIGFzc3VtcHRpb25zIGJlaGluZCB0aGUgc3RhdGlzdGljYWwgbW9kZWwgeW91J2xsIGZpdD8KCiAgICAxLiAgQXJlIHRob3NlIGFzc3VtcHRpb25zIHNhdGlzZmllZD8KCjYuICBGaXQgdGhlIG1vZGVsCgogICAgMS4gIEhvdyB3aWxsIHlvdSBhc3Nlc3Mgd2hldGhlciB0aGUgbW9kZWwgZml0cyB3ZWxsPwoKICAgIDIuICBDYW4geW91IGRldGVjdCBhbiBlZmZlY3Qgb2YgdGhlIHByZWRpY3RvcnM/CgogICAgMy4gIEhvdyBkbyB5b3UgbWVhc3VyZSB0aGUgZWZmZWN0PwoKNy4gIFdoYXQgZG8geW91IGNvbmNsdWRlIChpbmNsdWRpbmcgYW55IGNhdXRpb25zKQoKLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tLS0tCgojIyBBCgpAZGl4b25EaXNwcm9wb3J0aW9uYXRlSW1wb3J0YW5jZUxvbmcyMDE4IGZvY3VzZWQgb24gYXNzZW1ibGFnZXMgb2YgcmVwdGlsZXMgaW4gZm9yZXN0cyBhbmQgd29vZGxhbmRzIG9mIHNvdXRoZWFzdGVybiBBdXN0cmFsaWEsIHdpdGggYSBwYXJ0aWN1bGFyIGludGVyZXN0IGluIHJlbGF0aW9uc2hpcHMgYmV0d2VlbiB0aGVzZSBhc3NlbWJsYWdlcyBhbmQgdGhlIGZpcmUgaGlzdG9yeSBvZiB0aGUgbGFuZHNjYXBlLiBUaGV5IGlkZW50aWZpZWQgODEgc2l0ZXMgdGhhdCB2YXJpZWQgaW4gdGltZSBzaW5jZSBmaXJlLCBmcm9tIDYgbW9udGhzIHRvIFw+OTYgeS4gUmF0aGVyIHRoYW4gYSBjb250aW51dW0sIHRocmVlIGNhdGVnb3JpZXMgb2YgdGltZSBzaW5jZSBmaXJlIHdlcmUgdXNlZCwgMC41LTJ5LCA2LTEyeSwgYW5kIFw+OTZ5LiBTaXRlcyB3ZXJlIGFsc28gY2xhc3NpZmllZCBhY2NvcmRpbmcgdG8gaGFiaXRhdC4KClJlcHRpbGUgYXNzZW1ibGFnZXMgd2VyZSBzYW1wbGVkIHdpdGggYSByYW5nZSBvZiBtZXRob2RzLCBpbmNsdWRpbmcgdmlzdWFsIHN1cnZleXMgYW5kIGNhbWVyYSB0cmFwcywgdG8gZ2l2ZSBjb3VudHMgb2YgMjAgcmVwdGlsZXMsIHdob3NlIGFidW5kYW5jZSByYW5nZWQgZnJvbSAwLTEgdG8gMC0xMjYsIGRlcGVuZGluZyBvbiBzcGVjaWVzLgoKRGF0YSBhcmUgYXZhaWxhYmxlIGZyb20gW2RyeWFkXShodHRwczovL2RhdGFkcnlhZC5vcmcvc3Rhc2gvZGF0YXNldC9kb2k6MTAuNTA2MS9kcnlhZC5naDVtNDcxKSwgYXMgUmVwX2FidW5kLmNzdi4gWW91J2xsIHdhbnQgdG8gZm9jdXMgb24gdGhlIGNvbHVtbnMgd2l0aCByZXB0aWxlIG51bWJlcnMsIGFuZCB0aGUgb25lIHdpdGggdGhlIGZpcmUgaGlzdG9yeSAodHNmKS4gRm9yIGNvbnZlbmllbmNlLCB3ZSd2ZSBleHRyYWN0ZWQgdGhvc2UgZGF0YSBhcyBkaXhvbmJpb3RhLmNzdgoKKipTdGFydCBieSBoYXZpbmcgYSBxdWljayBsb29rIGF0IHRoZSBkYXRhLioqCgpgYGB7cn0KZGYgPC0gcmVhZF9jc3YoImRhdGEvZGl4b25iaW90YS5jc3YiKQpoZWFkKGRmLCAxMCkKYGBgCgpUaGUgZmlsZSBpcyBwcmV0dHkgc3RyYWlnaHRmb3J3YXJkLCB3aXRoIHRzZiBiZWluZyB0aGUgZmlyZSBjYXRlZ29yeSwgYW5kIGNvbHVtbnMgMy0yMiBlYWNoIHJlcHJlc2VudGluZyBvbmUgcmVwdGlsZSB0YXhvbi4KCiMjIE91dGxpbmUgaG93IHlvdSdkIGFzc2VzcyB3aGV0aGVyIHRoZSByZXB0aWxlIGFzc2VtYmxhZ2UgKHRoZSBjb21iaW5hdGlvbiBvZiBzcGVjaWVzIHByZXNlbnQgYW5kIHRoZWlyIGFidW5kYW5jZXMpIGlzIHJlbGF0ZWQgdG8gZmlyZSBoaXN0b3J5LCB1c2luZyBhIGRpc3NpbWlsYXJpdHktYmFzZWQgYXBwcm9hY2guCgpZb3Ugc2hvdWxkIHRoaW5rIGFib3V0IGhvdyB5b3UnbGwgZGVhbCB3aXRoIHRoZSBkYXRhOgoKLSAgIFdpbGwgeW91IG5lZWQgYSB0cmFuc2Zvcm1hdGlvbiwgc3RhbmRhcmRpemF0aW9uLCBldGMuPwoKLSAgIFdoYXQgYXJlIHRoZSBjb25zZXF1ZW5jZXMgb2YgYW55IGRlY2lzaW9ucyB5b3UgbWFrZSBmb3IgdGhlIGludGVycHJldGF0aW9uIG9mIHlvdXIgYW5hbHlzaXM/CgotICAgV2hhdCBtZWFzdXJlKHMpIG9mIGRpc3NpbWlsYXJpdHkgYXJlIGFwcHJvcHJpYXRlPwoKIyMjIFJ1biB0aGUgYW5hbHlzaXMgYW5kIHByb3ZpZGUgeW91ciBpbnRlcnByZXRhdGlvbgoKPiAjIyMgTURTIG9uIHJlcHRpbGUgYWJ1bmRhbmNlcwoKPiBTdGFydCB3aXRoIGRpc3NpbWlsYXJpdGllcwo+Cj4gU3R1ZGVudHMgc2hvdWxkIGhhdmUgc29tZSBkaXNjdXNzaW9uIG9mIHJhdyB2cyB0cmFuc2Zvcm1lZCAoZS5nLiA0dGggcm9vdCksIHN0YW5kYXJkaXphdGlvbi4KPgo+IEZvdXJ0aCByb290IGlzIGNvbW1vbi4gVXNpbmcgaXQgZG93bndlaWdodHMgdGhlIGluZmx1ZW5jZSBvZiBhYnVuZGFudCBzcGVjaWVzLCB3aGlsZSByYXcgYWxsb3dzIHRob3NlIHNwZWNpZXMgdG8gaGF2ZSBhIHN0cm9uZyBpbmZsdWVuY2UuIFRoZXJlIGFyZSBxdWVzdGlvbnMgYWJvdXQgd2hpY2ggb25lIG9mIHRoZXNlIHJlZmxlY3RzIGVjb2xvZ2ljYWwgZnVuY3Rpb24gYmV0dGVyCj4KPiBEYXRhIGFyZSBhbGwgYWJ1bmRhbmNlcwoKYGBge3J9CmRmMSA8LSBkZlssLSgxOjIpXSAgICNTdHJpcHBlZCBiYWNrIGZpbGUgd2l0aG91dCBsYWJlbHMsIGxlYXZpbmcgb25seSAKZGYxcyA8LSB3aXNjb25zaW4oZGYxKSAgICNDb21tb24gdG8gZG8gc29tZSBzdGFuZGFyZGlzYXRpb24KZGYxcy5iYyA8LSB2ZWdkaXN0KGRmMXMsJ2JyYXknKSAgICNXaXRoIHN0YW5kYXJkaXNhdGlvbiAod2lzY29uc2luKQpgYGAKCj4gVmlzdWFsaXplIHJlc3VsdHMKCj4gTmVlZCB0byBkZWNpZGUgd2hldGhlciB0byB1c2UgMiBvciAzIGRpbWVuc2lvbnMKCmBgYHtyfQojUnVuIDNkIGZpcnN0CmRmMXMubWRzMyA8LSBtZXRhTURTKGRmMXMuYmMsaz0zLGF1dG90cmFuc2Zvcm09RkFMU0UsdHJ5PTQwLHRyeW1heD04MCxtYXhpdD0yMDApCnN0cmVzc3Bsb3QoZGYxcy5tZHMzLCBtYWluPSJTaGVwYXJkIHBsb3QiKQpkZjFzLm1kczMKIyBOb3cgbG9vayBhdCAyZApkZjFzLm1kczIgPC0gbWV0YU1EUyhkZjFzLmJjLGs9MixhdXRvdHJhbnNmb3JtPUZBTFNFLHRyeT00MCx0cnltYXg9ODAsbWF4aXQ9MjAwKQpzdHJlc3NwbG90KGRmMXMubWRzMiwgbWFpbj0iU2hlcGFyZCBwbG90IikKZGYxcy5tZHMyCmBgYAoKPiBTdHJlc3MgdmFsdWVzIHN1Z2dlc3QgM2QgcGxvdCBpcyBhIG1vcmUgYWNjdXJhdGUgcmVwcmVzZW50YXRpb24KCj4gTmVlZCB0byBwbG90IGRpbSAxIHZzIGRpbSAyIGFuZCBkaW0gMSB2cyBkaW0gMyB3aXRoIHRzZiBncm91cHMgaWRlbnRpZmllZCBieSBzeW1ib2xzCj4KPiAqKkdlbmVyYXRlIGdyYXBocyoqCj4KPiBUaGlzIGlzIGEgdmFyaWFudCBvZiBjb2RlIHVzZWQgaW4gd29ya2VkIGV4YW1wbGUgZm9yIENoIDE2LiBJdCBkb2VzIHVzZSB0aGUgZ3JhcGggc2V0dGluZ3MgZGVmaW5lZCBpbiAqYXBwZWFyYW5jZS5SKiwgd2hpY2ggaXMgY2FsbGVkIGVhcmxpZXIuCgpgYGB7cn0KYTwtYXMuZGF0YS5mcmFtZShzY29yZXMoZGYxcy5tZHMzKSkKYTwtY2JpbmQoZGZbYygxOjIpXSxhKSAgICNBZGQgc2l0ZSBuYW1lcyAmIHN5bWJvbHMgZnJvbSBvcmlnaW5hbCBkYXRhIGZpbGUKcDE8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvcj10c2YsICkgKSsKICBnZW9tX3BvaW50KCkrCiAgbGFicyh5PSJNRFMyIiwgeD0iTURTMSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCh2YWx1ZXM9c3ltMywKICAgICAgICAgICAgICAgICAgICAgbmFtZT0iVFNGIiwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcDI8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMzLCBjb2xvciA9IHRzZiwgKSApKwogIGdlb21fcG9pbnQoKSsKICBsYWJzKHk9Ik1EUzMiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1zeW0zLAogICAgICAgICAgICAgICAgICAgICBuYW1lPSJUU0YiLAogICAgICAgICAgICAgICAgICAgICBndWlkZSA9CiAgICAgICAgICAgICAgICAgICAgICAgICBndWlkZV9sZWdlbmQobGFiZWwudGhlbWUgPSBlbGVtZW50X3RleHQoc2l6ZT02KSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGl0bGU9TlVMTCkKICAgICAgICAgICAgICAgICAgICAgKQoKcDEgKyBwMiArIHBsb3RfbGF5b3V0KGd1aWRlcz0nY29sbGVjdCcpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKYGBgCgojIyMgSG93IHdvdWxkIHlvdSBmaXQgYSBsaW5lYXIgbW9kZWwgdG8gYXNzZXNzIGZpcmUgZWZmZWN0cz8KCj4gRG8gcGVybWFub3ZhIGluY2x1ZGluZyB0c2YKCmBgYHtyfQpkZjFzLmFkbyA8LSBhZG9uaXMyKGRmMXMuYmN+dHNmLGRhdGE9ZGYscGVybXV0YXRpb25zPTk5OSkKcHJpbnQoZGYxcy5hZG8pCmBgYAoKPiBQZXJtYW5vdmEgc3VnZ2VzdHMgc3Ryb25nIHNpZ25hbCBmcm9tIHRzZgoKPiBTSU1QRVIgYW5hbHlzaXMgZm9yIHRzZiBncm91cCBkaWZmZXJlbmNlcwoKYGBge3J9CmRmLnNpbSA8LSBzaW1wZXIoZGYxc1ssLSgxOjIpXSwgZGYkdHNmLCBwZXJtdXRhdGlvbnM9MTAwMCkKc3VtbWFyeShkZi5zaW0pCmBgYAoKIyMgU3VwcG9zZSB5b3UgdGhpbmsgb2YgYSByZXB0aWxlIGFzc2VtYmxhZ2UgYXMgc2ltcGx5IHRoZSBzcGVjaWVzIHRoYXQgYXJlIHByZXNlbnQsIHJlZ2FyZGxlc3Mgb2YgaG93IGNvbW1vbiB0aGV5IGFyZS4gV2hhdCB3b3VsZCB5b3UgY29uY2x1ZGUgYWJvdXQgdGhlIGluZnVlbmNlIG9mIHRpbWUgc2luY2UgZmlyZSBub3c/Cgo+IFJlcGVhdCB0aGUgcHJlY2VkaW5nIGFuYWx5c2lzIGFmdGVyIGNyZWF0aW5nIGEgbmV3IGRhdGEgZmlsZSB0aGF0IGlzIHByZXNlbmNlLWFic2VuY2Ugb25seS4KCmBgYHtyfQpkZjEuamFjIDwtIHZlZ2Rpc3QoZGYxLGJpbmFyeT1UUlVFLG1ldGhvZD0namFjY2FyZCcpCiNSdW4gM2QgZmlyc3QKZGYxcy5tZHMzIDwtIG1ldGFNRFMoZGYxLmphYyxrPTMsYXV0b3RyYW5zZm9ybT1GQUxTRSx0cnk9NDAsdHJ5bWF4PTgwLG1heGl0PTIwMCkKc3RyZXNzcGxvdChkZjFzLm1kczMsIG1haW49IlNoZXBhcmQgcGxvdCIpCmRmMXMubWRzMwojIE5vdyBsb29rIGF0IDJkCmRmMXMubWRzMiA8LSBtZXRhTURTKGRmMS5qYWMsaz0yLGF1dG90cmFuc2Zvcm09RkFMU0UsdHJ5PTQwLHRyeW1heD04MCxtYXhpdD0yMDApCnN0cmVzc3Bsb3QoZGYxcy5tZHMyLCBtYWluPSJTaGVwYXJkIHBsb3QiKQpkZjFzLm1kczIKYGBgCgo+IE1hcmdpbmFsIGRlY2lzaW9uIGhlcmUgLSBzdHJlc3MganVzdCBvdmVyIDAuMiBmb3Igaz0yLCBzbyBtaWdodCBiZSBPSy4gU3RyZXNzIGdvb2QgZm9yIGs9MywgYnV0IDItZCBNRFMgdXN1YWxseSBiZXR0ZXIgZm9yIHJlcG9ydGluZyBvciBzaG93aW5nIHRvIGF1ZGllbmNlCgpgYGB7cn0KYTwtYXMuZGF0YS5mcmFtZShzY29yZXMoZGYxcy5tZHMzKSkKYTwtY2JpbmQoZGZbYygxOjIpXSxhKSAgICNBZGQgc2l0ZSBuYW1lcyAmIHN5bWJvbHMgZnJvbSBvcmlnaW5hbCBkYXRhIGZpbGUKcDE8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvcj10c2YsICkgKSsKICBnZW9tX3BvaW50KCkrCiAgbGFicyh5PSJNRFMyIiwgeD0iTURTMSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCh2YWx1ZXM9c3ltMywKICAgICAgICAgICAgICAgICAgICAgbmFtZT0iVFNGIiwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcDI8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMzLCBjb2xvciA9IHRzZiwgKSApKwogIGdlb21fcG9pbnQoKSsKICBsYWJzKHk9Ik1EUzMiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1zeW0zLAogICAgICAgICAgICAgICAgICAgICBuYW1lPSJUU0YiLAogICAgICAgICAgICAgICAgICAgICBndWlkZSA9CiAgICAgICAgICAgICAgICAgICAgICAgICBndWlkZV9sZWdlbmQobGFiZWwudGhlbWUgPSBlbGVtZW50X3RleHQoc2l6ZT02KSwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGl0bGU9TlVMTCkKICAgICAgICAgICAgICAgICAgICAgKQoKcDEgKyBwMiArIHBsb3RfbGF5b3V0KGd1aWRlcz0nY29sbGVjdCcpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKYGBgCgo+IFVuYnVybmVkIHNpdGVzIHNlcGFyYXRlIHZlcnkgY2xlYXJseSBvbiBNRFMxLiBQYXR0ZXJuIGNsZWFyZXIgdGhhbiB3aXRoIGFidW5kYW5jZXMgaW5jbHVkZWQKCj4gRG8gcGVybWFub3ZhIGluY2x1ZGluZyB0c2YKCmBgYHtyfQpkZjEuYWRvIDwtIGFkb25pczIoZGYxLmphY350c2YsZGF0YT1kZixwZXJtdXRhdGlvbnM9OTk5KQpwcmludChkZjEuYWRvKQpgYGAKCj4gQ2xlYXIgZmlyZSBzaWduYWwKCiMjIEV4dGVuc2lvbiBhY3Rpdml0eQoKRGl4b24gYW5kIGhlciBjb2xsZWFndWVzIGFsc28gcmVjb3JkZWQgYSByYW5nZSBvZiBoYWJpdGF0IHZhcmlhYmxlcywgd2hpY2ggd2UndmUgY29sbGVjdGVkIGZvciB5b3UgaW4gZGl4b25lbnYuY3N2CgpgYGB7cn0KZGl4b25lbnYgPC0gcmVhZF9jc3YoImRhdGEvZGl4b25lbnYuY3N2IikKaGVhZChkaXhvbmVudiwgMTApCmBgYAoKVGhleSByZWNvcmRlZCBDb2Fyc2UgV29vZHkgRGVicmlzLCB3aGljaCB3YXMgbG9uZy10cmFuc2Zvcm1lZCwgbGl0dGVyIGNvdmVyIChsaXRjb3YpLCAlIGNvdmVyIG9mIGdyb3VuZGNvdmVyIChncmNvdiksIGFuZCAlIGNvdmVyIG9mIHNocnVicyAoc2hyY292KSwgYW5kIHJvY2sgY292ZXIuIEluIHRoZSBvcmlnaW5hbCBwYXBlciwgRGl4b24gZXQgYWwuIHdlcmUgY29tZm9ydGFibGUgdGhhdCB0aGVzZSBwcmVkaWN0b3JzIHdlcmVuJ3QgY29ycmVsYXRlZC4KCiMjIyBIb3cgaXMgcmVwdGlsZSBhc3NlbWJsYWdlIHJlbGF0ZWQgdG8gdGltZSBzaW5jZSBmaXJlIGFuZCB0aGVzZSBmaXZlIGhhYml0YXQgdmFyaWFibGVzCgo+IENoZWNrIGNvbnRpbnVvdXMgcHJlZGljdG9ycwoKYGBge3J9CmJveHBsb3QoZGl4b25lbnYkbGN3ZCkKYm94cGxvdChkaXhvbmVudiRsaXRjb3YpCmJveHBsb3QoZGl4b25lbnYkZ3Jjb3YpCmJveHBsb3QoZGl4b25lbnYkc2hyY292KQpib3hwbG90KGRpeG9uZW52JHJvY2tzKQojdHJhbnNmb3JtIHJvY2tzIHRvIGxvZysxCmRpeG9uZW52JGxyb2NrcyA8LSBsb2cxMChkaXhvbmVudiRyb2NrcysxKQojIGNoZWNrIGhvbW9nIG9mIGRpc3BlcnNpb25zIG9uIGRvdWJsZSBzdGFuZGFyZGl6ZWQgYWJ1bmRhbmNlcwpkaXhvbmJpb3RhMXMuZGlzcCA8LSBiZXRhZGlzcGVyKGRpeG9uYmlvdGExcy5iYyxkaXhvbmJpb3RhJHRzZikKYW5vdmEoZGl4b25iaW90YTFzLmRpc3ApCmBgYAoKPiBkbyBwZXJtYW5vdmEgaW5jbHVkaW5nIHRzZiBwbHVzIGNvbnRpbnVvdXMgcHJlZGljdG9ycwoKYGBge3J9CmRpeG9uYmlvdGExcy5hZG8gPC0gYWRvbmlzMihkaXhvbmJpb3RhMXMuYmN+dHNmK2xjd2QrbGl0Y292K2dyY292K3NocmNvditscm9ja3MsZGF0YT1kaXhvbmVudixwZXJtdXRhdGlvbnM9OTk5KQpwcmludChkaXhvbmJpb3RhMXMuYWRvKQojUnVuIGFuYWx5c2lzIGZvciBwcmVzZW5jZS1hYnNlbmNlIGFzIHdlbGwKZGl4b25iaW90YTEuamFjLmFkbyA8LSBhZG9uaXMyKGRpeG9uYmlvdGExLmphY350c2YrbGN3ZCtsaXRjb3YrZ3Jjb3Yrc2hyY292K2xyb2NrcyxkYXRhPWRpeG9uZW52LHBlcm11dGF0aW9ucz05OTkpCnByaW50KGRpeG9uYmlvdGExLmphYy5hZG8pCgpgYGAKCj4gTWFpbiBlZmZlY3QgaXMgdHNmOyBvbmx5IHNvbWUgcXVlc3Rpb24gYWJvdXQgQ1dEIHdoZW4gd2UgbG9vayBhdCBzcGVjaWVzIGNvbXBvc2l0aW9uIGFuZCBhYnVuZGFuY2UsIHdoaWxlIGZvciBwcmVzZW5jZS1hYnNlbmNlIGRhdGEsIHRoZXJlIGFyZSBhbHNvIGVmZmVjdHMgb2YgQ1dEIGFuZCBncm91bmRjb3ZlcgoKIyMjIG12YWJ1bmQKCmBgYHtyIGVycm9yPVRSVUV9CmRpeG9uYmlvdGFtdiA8LSBtdmFidW5kKGRpeG9uYmlvdGFbLC0oMToyKV0pCmRpeG9uYmlvdGFtdi5tdiA8LSBtYW55Z2xtKGRpeG9uYmlvdGFtdn5kaXhvbmVudiR0c2YrZGl4b25lbnYkbGN3ZCtkaXhvbmVudiRsaXRjb3YrZGl4b25lbnYkZ3Jjb3YrZGl4b25lbnYkc2hyY292K2RpeG9uZW52JGxyb2NrcyxmYW1pbHk9InBvaXNzb24iKQpwbG90KGRpeG9uYmlvdGFtdi5tdikKIyBzdGlsbCBoZXRlcm8gdmFycyBzbyB1c2UgLXZlIGJpbm9taWFsCmRpeG9uYmlvdGFtdjEubXYgPC0gbWFueWdsbShkaXhvbmJpb3RhbXZ+ZGl4b25lbnYkdHNmK2RpeG9uZW52JGxjd2QrZGl4b25lbnYkbGl0Y292K2RpeG9uZW52JGdyY292K2RpeG9uZW52JHNocmNvditkaXhvbmVudiRscm9ja3MsZmFtaWx5PSJuZWdhdGl2ZV9iaW5vbWlhbCIpCnBsb3QoZGl4b25iaW90YW12MS5tdikKYW5vdmEoZGl4b25iaW90YW12MS5tdikKYGBgCgojIyBCIFNpbXBsZSBNRFMgJiBQZXJtYW5vdmEKCkxldCdzIHJldHVybiB0byB0aGUgQGh1dHRvVXJiYW5PcGVuU3BhY2VzMjAyMSBleGFtcGxlIGZyb20gdGhlIENoYXB0ZXIgMTUgZXhlcmNpc2VzLiBUaGVyZSwgd2UgdXNlZCBzaW1pbGFyaXR5LWJhc2VkIGFuYWx5c2VzIHRvIGV4cGxvcmUgdGhlIHJlbGF0aW9uc2hpcCBiZXR3ZWVuIGZyb2cgYXNzZW1ibGFnZXMgaW4gcG9uZHMgYW5kIHRoZSBkZWdyZWUgb2YgdXJiYW5pemF0aW9uIHN1cnJvdW5kaW5nLiBUaGlzIHNlZW1zIGEgcXVlc3Rpb24gdGhhdCBjb3VsZCBqdXN0IGFzIGVhc2lseSBiZSBleGFtaW5lZCB1c2luZyBkaXN0YW5jZXMgb3IgZGlzc2ltaWxhcml0aWVzLgoKUmV0dXJuIHRvIHRoZSBkYXRhIGFuZCBhc3Nlc3Mgd2hldGhlciBmcm9nIGFzc2VtYmxhZ2VzICh1c2luZyB0aGUgc3RhbmRhcmRpemVkIGFidW5kYW5jZSBzY2FsZSBvciBwcmVzZW5jZS1hYnNlbmNlKSBkaWZmZXIgd2l0aCB1cmJhbml6YXRpb24uCgojIyMjIERpZCB5b3VyIGNvbmNsdXNpb25zIGRpZmZlciBmcm9tIHRob3NlIG9idGFpbmVkIHVzaW5nIFJEQT8KCmBgYHtyfQpkZiA8LSByZWFkX2NzdigiZGF0YS9odXR0b2FtcGguY3N2IikKaGVhZChkZiwxMCkKI3BvbmRzIDkgYW5kIDQwIGhhZCBubyBmcm9ncy4gUmVtb3ZlIGZvciBkaXN0YW5jZSBhbmFseXNpcwpkZiA8LSBkZiAlPiUgCiAgZmlsdGVyKFNpdGUhPSA5ICYgU2l0ZSAhPTQwKQpkZgpgYGAKClRoZSBleHBsYW5hdGlvbiBvZiB0aGUgdmFyaWFibGVzIGlzIGluIHRoZSBwcmV2aW91cyBjaGFwdGVyJ3MgZXhlcmNpc2VzLgoKYGBge3J9CmRmMSA8LSBkZlssLSgxOjIpXSAgICNTdHJpcHBlZCBiYWNrIGZpbGUgd2l0aG91dCBsYWJlbHMsIGxlYXZpbmcgb25seSAKZGYxLmJjIDwtIHZlZ2Rpc3QoZGYxLCdicmF5JykKYGBgCgpgYGB7cn0KI1J1biAyZCBmaXJzdApkZjEubWRzMiA8LSBtZXRhTURTKGRmMS5iYyxrPTIsYXV0b3RyYW5zZm9ybT1GQUxTRSx0cnk9NDAsdHJ5bWF4PTgwLG1heGl0PTIwMCkKc3RyZXNzcGxvdChkZjEubWRzMiwgbWFpbj0iU2hlcGFyZCBwbG90IikKZGYxLm1kczIKIyBOb3cgbG9vayBhdCAzZApkZi5tZHMzIDwtIG1ldGFNRFMoZGYxLmJjLGs9MyxhdXRvdHJhbnNmb3JtPUZBTFNFLHRyeT00MCx0cnltYXg9ODAsbWF4aXQ9MjAwKQpzdHJlc3NwbG90KGRmLm1kczMsIG1haW49IlNoZXBhcmQgcGxvdCIpCmRmLm1kczMKYGBgCgo+IDItZCBpcyBhY2NlcHRhYmxlCgpgYGB7cn0KYTwtYXMuZGF0YS5mcmFtZShzY29yZXMoZGYxLm1kczIpKQphPC1jYmluZChkZltjKDE6MildLGEpICAgI0FkZCBzaXRlIG5hbWVzICYgc3ltYm9scyBmcm9tIG9yaWdpbmFsIGRhdGEgZmlsZQpnZ3Bsb3QoZGF0YT1hLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3I9dHlwZSwgKSApKwogIGdlb21fcG9pbnQoKSsKICBsYWJzKHk9Ik1EUzIiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1zeW0zLAogICAgICAgICAgICAgICAgICAgICBuYW1lPSJVcmJhbml6YXRpb24iLAogICAgICAgICAgICAgICAgICAgICBndWlkZSA9CiAgICAgICAgICAgICAgICAgICAgICAgICBndWlkZV9sZWdlbmQobGFiZWwudGhlbWUgPSBlbGVtZW50X3RleHQoc2l6ZT02KSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGl0bGU9TlVMTCkKICAgICAgICAgICAgICAgICAgICAgKSsKICB0aGVtZV9xaygpICsgc2NhbGVfY29sb3JfdWNoaWNhZ28oKQpgYGAKCmBgYHtyfQpkZjEuYWRvIDwtIGFkb25pczIoZGYxLmJjfnR5cGUsZGF0YT1kZixwZXJtdXRhdGlvbnM9OTk5KQpwcmludChkZjEuYWRvKQpgYGAKCj4gUGVybWFub3ZhIGRldGVjdHMgYSBkaWZmZXJlbmNlOyBNRFMgcGxvdCBzdWdnZXN0cyBsb3cgdXJiYW5pemF0aW9uIG1heSBiZSBkaWZmZXJlbnQgZnJvbSB0aGUgb3RoZXIgdHdvLCB0aG91Z2ggd29ydGggbG9va2luZyBhdCBkaXNwZXJzaW9uIGFzIHdlbGwsIGFzIHNwcmVhZCBvZiB0aGlzIHR5cGUgaXMgbGVzcyB0aGFuIHRoZSBvdGhlciB0d28gb24gcGxvdAoKPiBSZXBlYXQgYW5hbHlzaXMgd2l0aCBwcmVzZW5jZS1hYnNlbmNlCgpkZjEuamFjIFw8LSB2ZWdkaXN0KGRmMSxiaW5hcnk9VFJVRSxtZXRob2Q9J2phY2NhcmQnKQoKYGBge3J9CiNSdW4gMmQgZmlyc3QKZGYxLmphYyA8LSB2ZWdkaXN0KGRmMSxiaW5hcnk9VFJVRSxtZXRob2Q9J2phY2NhcmQnKQpkZjEubWRzMiA8LSBtZXRhTURTKGRmMS5qYWMsaz0yLGF1dG90cmFuc2Zvcm09RkFMU0UsdHJ5PTQwLHRyeW1heD04MCxtYXhpdD0yMDApCnN0cmVzc3Bsb3QoZGYxLm1kczIsIG1haW49IlNoZXBhcmQgcGxvdCIpCmRmMS5tZHMyCiMgTm93IGxvb2sgYXQgM2QKZGYubWRzMyA8LSBtZXRhTURTKGRmMS5qYWMsaz0zLGF1dG90cmFuc2Zvcm09RkFMU0UsdHJ5PTQwLHRyeW1heD04MCxtYXhpdD0yMDApCnN0cmVzc3Bsb3QoZGYubWRzMywgbWFpbj0iU2hlcGFyZCBwbG90IikKZGYubWRzMwpgYGAKCj4gQWdhaW4sIHN0cmVzcyBPSyB3aXRoIDJELCBnb29kIHdpdGggM0QuIExvb2sgYXQgMkQgZmlyc3QKCmBgYHtyfQphPC1hcy5kYXRhLmZyYW1lKHNjb3JlcyhkZjEubWRzMikpCmE8LWNiaW5kKGRmW2MoMToyKV0sYSkgICAjQWRkIHNpdGUgbmFtZXMgJiBzeW1ib2xzIGZyb20gb3JpZ2luYWwgZGF0YSBmaWxlCmdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvcj10eXBlLCApICkrCiAgZ2VvbV9wb2ludCgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiKSsKICBzY2FsZV9zaGFwZV9tYW51YWwodmFsdWVzPXN5bTMsCiAgICAgICAgICAgICAgICAgICAgIG5hbWU9IlVyYmFuaXphdGlvbiIsCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApKwogIHRoZW1lX3FrKCkgKyBzY2FsZV9jb2xvcl91Y2hpY2FnbygpCmBgYAoKYGBge3J9CmRmMS5hZG8gPC0gYWRvbmlzMihkZjEuamFjfnR5cGUsZGF0YT1kZixwZXJtdXRhdGlvbnM9OTk5KQpwcmludChkZjEuYWRvKQpgYGAKCj4gTm8gc2VwYXJhdGlvbiBzZWVuIG9uIFByZXNlbmNlLUFic2VuY2UKCiMjIEMKCkBncmlmZmVuU2lnbmlmaWNhbnRFZmZlY3RISVYyMDE5IGV4YW1pbmVkIHRoZSBvcmFsIG1pY3JvYmlvbWUgb2YgaHVtYW5zLCB3aXRoIGEgZm9jdXMgb24gdGhlIGVmZmVjdHMgb2YgSElWIGFuZCBBbnRpUmV0cm92aXJhbCB0aGVyYXB5IChBUlQpIG9uIHRoaXMgbWljcm9iaW9tZS4gVGhleSBkZXNjcmliZWQgdGhlIG1pY3JvYmlvbWUgb2YgMzQxIHBhdGllbnRzLCB3aG8gZmVsbCBpbnRvIG9uZSBvZiB0d28gY2F0ZWdvcmllczogSElWXi1eLCBhbmQgSElWKyB3aXRoIEFSVC4gVGhleSBhbHNvIG1hdGNoZWQgdGhlaXIgc2FtcGxlcyBhcyBmYXIgYXMgcG9zc2libGUgZm9yIHNleCwgYWxvbmcgd2l0aCBvdGhlciBjb25kaXRpb25zIHRoYXQgY2FuIGluZmx1ZW5jZSBvcmFsIG1pY3JvYmlvbWVzICgqQ2FuZGlkYSogaW5mZWN0aW9uLCBjdXJyZW50IHNtb2tpbmcpLiBXZSdsbCB1c2UgdGhlaXIgcmVjb3JkcyBmb3IgdGhlc2Ugb3RoZXIgY2F0ZWdvcmllcyBhcyB3ZWxsLiBUaGVyZSB3ZXJlIG1vcmUgSElWKyB0aGFuIC0gc3ViamVjdHMsIGJ1dCB3aXRoaW4gdGhlc2UgZ3JvdXBzLCBhcHByb3hpbWF0ZWx5IGVxdWFsIG51bWJlcnMgb2YgdHdvIHNleGVzLCB0d28gKkNhbmRpZGEqIGluZmVjdGlvbiBzdGF0dXMsIGFuZCB0d28gc21va2luZyBjYXRlZ29yaWVzLgoKVGhlIGJhY3RlcmlvbWUgd2FzIHJlY29yZGVkIHNlcGFyYXRlbHkgdXNpbmcgMTZTIFJOQSBzZXF1ZW5jaW5nLCB3aGljaCBvdmVyZWQganVzdCBvdmVyIDYwMCB0YXhhLgoKVGhlIGRhdGEgYXJlIGF2YWlsYWJsZSBmcm9tIHRoZWlyIHBhcGVyLCBhcyB0d28gRXhjZWwgZmlsZXMgaW4gdGhlIHN1cHBsZW1lbnRhcnkgaW5mb3JtYXRpb24uIFRoZSBtZXRhZGF0YSBnaXZlcyBISVYgc3RhdHVzLCBhbmQgYSByYWZ0IG9mIGRlbW9ncmFwaGljIGluZm9ybWF0aW9uLiBXZSdsbCBqdXN0IHdvcmsgd2l0aCBISVYsIHNleCwgKkNhbmRpZGEqLCBhbmQgY3VycmVudCBzbW9raW5nLiBBIHNlY29uZCBzaGVldCBoYXMgdGhlIGJhY3RlcmlvbWUuIFRoZSBjb2RlIGNodW5rIGJlbG93IHJlYWRzIHRoaXMgZmlsZSAoZnJvbSBhIERvd25sb2FkcyBmb2xkZXIgLSB5b3UnbGwgbmVlZCB0byBtb2RpZnkgdGhlIGZpbGUgbG9jYXRpb24pLiBJdCBhbHNvIHNjcmVlbnMgdGhlIGZpbGUgZm9yIGFueSBiYWN0ZXJpYWwgdGF4YSB0aGF0IHdlcmUgbm90IHByZXNlbnQgaW4gdGhpcyBwYXJ0aWN1bGFyIHN0dWR5LCB3aGljaCByZWR1Y2VzIHRoZSBiYWN0ZXJpYWwgdGF4YSB0byA1OTkuCgpgYGB7ciBkZn0KbGlicmFyeShyZWFkeGwpCiNSZWFkIGluIG1ldGFkYXRhLiBMb3RzIG9mIGluZm9ybWF0aW9uIHdlIGRvbid0IG5lZWQsIHNvIGEgY291cGxlIG9mIGl0ZXJhdGlvbnMgdG8gZ2V0IGRvd24gdG8gYSBmZXcgY29sdW1ucyB3ZSB3YW50IC0gSElWIHN0YXR1cywgY2FuZGlkYSwgY3VycmVudCBzbW9raW5nLCBnZW5kZXIKZGYgPC0gcmVhZF9leGNlbCgifi9Eb3dubG9hZHMvNDE1OThfMjAxOV81NTcwM19NT0VTTTNfRVNNLnhsc3giLCAKICAgICByYW5nZSA9ICJBMTpTMzQyIiwgbmEgPSAiTkEiKQpkZjwtc2VsZWN0KGRmLCBISVYsIGNhbmRpZGEsIGdlbmRlciwgc21va2luZ19jdXJyZW50KQpoZWFkKGRmKQojZGYxIGlzIHRoZSBiYWN0ZXJpYWwgZGF0YQpkZjEgPC0gcmVhZF9leGNlbCgifi9Eb3dubG9hZHMvNDE1OThfMjAxOV81NTcwM19NT0VTTTJfRVNNLnhsc3giKQpoZWFkKGRmMSwxMCkKZGYxIDwtIGRmMSAlPiUKICBzZWxlY3Qod2hlcmUoIH4gaXMubnVtZXJpYygueCkgJiYgc3VtKC54KSAhPSAwKSkgICNEcm9wcyBhbnkgYmFjdGVyaWFsIHRheG9uIG5vdCByZWNvcmRlZCBpbiB0aGlzIHN0dWR5IChpLmUuIHdoZXJlIGl0J3MgYWxsIHplcm9lcykKZGYxIDwtIGRmMVssLTFdICNyZW1vdmUgY29sMSwgd2hpY2ggaXMgc2FtcGxlIElECmBgYAoKIyMjIFVzZSBhIGRpc3NpbWlsYXJpdHkgYmFzZWQgYXBwcm9hY2ggdG8gYXNzZXNzIHRoZSBjb21iaW5lZCBlZmZlY3RzIG9mIEhJVi1BUlQsIHNleCwgKkNhbmRpZGEqIGluZmVjdGlvbiwgYW5kIGN1cnJlbnQgc21va2luZyBvbiB0aGUgYmFjdGVyaW9tZQoKKipEbyB0aGVzZSBmYWN0b3JzIGFjdCBpbmRlcGVuZGVudGx5IG9uIHRoZSBiYWN0ZXJpb21lPyoqCgoqKldoaWNoIGVmZmVjdHMgYXJlIGxhcmdlc3QgKHVzZSBzb21lIGdyYXBoaWNhbCBtZXRob2RzKT8qKgoKKipGaXJzdCBzdGVwczoqKiBUaGluayBhYm91dCBzdGFuZGFyZGl6YXRpb24gYW5kIHdoaWNoIGRpc3RhbmNlIG1lYXN1cmUgeW91J2xsIHVzZS4gQnJheS1DdXJ0aXMgc2VlbXMgY29tbW9uIGZvciB0aGVzZSBraW5kIG9mIGRhdGEsIGFuZCBhbmQgY291bnRzIG9mIGRpZmZlcmVudCB0YXhhIHZhcnkgd2lkZWx5LgoKYGBge3IgQkN9CmRmMXMgPC0gd2lzY29uc2luKGRmMSkgICAjQ29tbW9uIHRvIGRvIHNvbWUgc3RhbmRhcmRpc2F0aW9uLCBhbmQgY291bnRzIG9mIGRpZmZlcmVudCB0YXhhIHZhcnkgd2lkZWx5CmRmMXMuYmMgPC0gdmVnZGlzdChkZjFzLCdicmF5JykgICAjV2l0aCBzdGFuZGFyZGlzYXRpb24gKHdpc2NvbnNpbikKYGBgCgo+UnVuIGZhY3RvcmlhbCBtb2RlbCB1c2luZyBwZXJtYW5vdmEsIGJhc2VkIG9uIEItQwoKYGBge3J9CmRmMS5hZG8gPC0gYWRvbmlzMihkZjFzLmJjfkhJVipjYW5kaWRhKmdlbmRlcipzbW9raW5nX2N1cnJlbnQsZGF0YT1kZixwZXJtdXRhdGlvbnM9OTk5KQpwcmludChkZjEuYWRvKQpgYGAKCj5Nb2RlbCBzaG93cyBtYWluIGVmZmVjdHMsIGJ1dCBubyBpbnRlcmFjdGlvbnMgaW1wb3J0YW50IChvciBhdCBsZWFzdCwgZGV0ZWN0ZWQpCgoqKkNvdWxkIHlvdSBmaXQgYSBzaW1wbGVyIG1vZGVsIHRvIHRoZSBkYXRhPyoqCgo+UnVuIHNpbXBsZXIgbW9kZWwKCmBgYHtyfQpkZjIuYWRvIDwtIGFkb25pczIoZGYxcy5iY35ISVYrY2FuZGlkYStnZW5kZXIrc21va2luZ19jdXJyZW50LGRhdGE9ZGYscGVybXV0YXRpb25zPTk5OSkKcHJpbnQoZGYyLmFkbykKYGBgCgo+Tm8gY2hhbmdlOyBub3QgdG9vIHN1cnByaXNpbmcsIGFzIGFsbCBmb3VyIG1haW4gZWZmZWN0cyB3ZXJlIGRldGVjdGVkIGVhc2lseS4KCj5Ob3cgZ28gdG8gZGF0YSB2aXN1YWxpemF0aW9uCgpgYGB7cn0KI1J1biAzZCBmaXJzdApkZjFzLm1kczMgPC0gbWV0YU1EUyhkZjFzLmJjLGs9MyxhdXRvdHJhbnNmb3JtPUZBTFNFLHRyeT00MCx0cnltYXg9ODAsbWF4aXQ9MjAwKQpzdHJlc3NwbG90KGRmMXMubWRzMywgbWFpbj0iU2hlcGFyZCBwbG90IikKZGYxcy5tZHMzCiMgTm93IGxvb2sgYXQgMmQKZGYxcy5tZHMyIDwtIG1ldGFNRFMoZGYxcy5iYyxrPTIsYXV0b3RyYW5zZm9ybT1GQUxTRSx0cnk9NDAsdHJ5bWF4PTgwLG1heGl0PTIwMCkKc3RyZXNzcGxvdChkZjFzLm1kczIsIG1haW49IlNoZXBhcmQgcGxvdCIpCmRmMXMubWRzMgpgYGAKCgo+M0QgbW9kZWwgZml0cyBiZXR0ZXIgdGhhbiAyRC4gU3RyZXNzIGZvciAyRCBhYm92ZSAwLjIKCj4gU3RhcnQgd2l0aCBwbG90cyB3aGVyZSBzeW1ib2wgY29sb3VyIGluZGljYXRlcyBsZXZlbHMgb2Ygb25lIG9mIHRoZSBmYWN0b3JzCgpgYGB7cn0KYTwtYXMuZGF0YS5mcmFtZShzY29yZXMoZGYxcy5tZHMzKSkKYTwtY2JpbmQoZGZbYygxOjQpXSxhKSAgICNBZGQgc2l0ZSBuYW1lcyAmIHN5bWJvbHMgZnJvbSBvcmlnaW5hbCBkYXRhIGZpbGUKcDE8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvcj1ISVYpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiKSsKICBzY2FsZV9zaGFwZV9tYW51YWwodmFsdWVzPWMoMTYsMTcpLAogICAgICAgICAgICAgICAgICAgICBndWlkZSA9CiAgICAgICAgICAgICAgICAgICAgICAgICBndWlkZV9sZWdlbmQobGFiZWwudGhlbWUgPSBlbGVtZW50X3RleHQoc2l6ZT02KSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGl0bGU9TlVMTCkKICAgICAgICAgICAgICAgICAgICAgKQpwMjwtZ2dwbG90KGRhdGE9YSwgYWVzKHg9Tk1EUzEsIHk9Tk1EUzMsIGNvbG9yID0gSElWKSApKwogIGdlb21fcG9pbnQoKSsKICBzdGF0X2VsbGlwc2UoKSsKICBsYWJzKHk9Ik1EUzMiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1jKDE2LDE3KSwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKCnAxICsgcDIgKyBwbG90X2xheW91dChndWlkZXM9J2NvbGxlY3QnKSArIHBsb3RfYW5ub3RhdGlvbih0aXRsZT0iSElWIHN0YXR1cyIpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKCnAxPC1nZ3Bsb3QoZGF0YT1hLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3I9Y2FuZGlkYSkgKSsKICBnZW9tX3BvaW50KCkrCiAgc3RhdF9lbGxpcHNlKCkrCiAgbGFicyh5PSJNRFMyIiwgeD0iTURTMSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCh2YWx1ZXM9YygxNiwxNyksCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApCnAyPC1nZ3Bsb3QoZGF0YT1hLCBhZXMoeD1OTURTMSwgeT1OTURTMywgY29sb3IgPSBjYW5kaWRhKSApKwogIGdlb21fcG9pbnQoKSsKICBzdGF0X2VsbGlwc2UoKSsKICBsYWJzKHk9Ik1EUzMiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1jKDE2LDE3KSwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKCnAxICsgcDIgKyBwbG90X2xheW91dChndWlkZXM9J2NvbGxlY3QnKSArIHBsb3RfYW5ub3RhdGlvbih0aXRsZT0iQ2FuZGlkYSIpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKCnAxPC1nZ3Bsb3QoZGF0YT1hLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3I9Z2VuZGVyKSApKwogIGdlb21fcG9pbnQoKSsKICBzdGF0X2VsbGlwc2UoKSsKICBsYWJzKHk9Ik1EUzIiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1jKDE2LDE3KSwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcDI8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMzLCBjb2xvciA9IGdlbmRlcikgKSsKICBnZW9tX3BvaW50KCkrCiAgc3RhdF9lbGxpcHNlKCkrCiAgbGFicyh5PSJNRFMzIiwgeD0iTURTMSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCh2YWx1ZXM9YygxNiwxNyksCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApCgpwMSArIHAyICsgcGxvdF9sYXlvdXQoZ3VpZGVzPSdjb2xsZWN0JykgKyBwbG90X2Fubm90YXRpb24odGl0bGU9IlNleCIpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKCnAxPC1nZ3Bsb3QoZGF0YT1hLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3I9c21va2luZ19jdXJyZW50KSApKwogIGdlb21fcG9pbnQoKSsKICBzdGF0X2VsbGlwc2UoKSsKICBsYWJzKHk9Ik1EUzIiLCB4PSJNRFMxIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKHZhbHVlcz1jKDE2LDE3KSwKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcDI8LWdncGxvdChkYXRhPWEsIGFlcyh4PU5NRFMxLCB5PU5NRFMzLCBjb2xvciA9IHNtb2tpbmdfY3VycmVudCkgKSsKICBnZW9tX3BvaW50KCkrCiAgc3RhdF9lbGxpcHNlKCkrCiAgbGFicyh5PSJNRFMzIiwgeD0iTURTMSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCh2YWx1ZXM9YygxNiwxNyksCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApCgpwMSArIHAyICsgcGxvdF9sYXlvdXQoZ3VpZGVzPSdjb2xsZWN0JykgKyBwbG90X2Fubm90YXRpb24odGl0bGU9IkN1cnJlbnQgc21va2luZyIpICYgdGhlbWVfcWsoKSAmIHNjYWxlX2NvbG9yX3VjaGljYWdvKCkKCmBgYAoKPlRyeSBzb21lIHZpc3VhbGl6YXRpb25zIG9mIHBhaXJzIG9mIGZhY3RvcnMgdXNpbmcgTURTMSAmIE1EUzIgKHdoaWxlIDNEIGlzIHByZWZlcnJlZCwgc3RyZXNzIG9mIDNEIGlzIGFyb3VuZCAuMTUsIHZzIDAuMjIsIHNvIG5vdCBodWdlIGRpZmZlcmVuY2UpCj5TZXBhcmF0ZSBISVYrLy0sIHRoZW4gdXNlIHN5bWJvbCBjb2xvdXIgdG8gc2VwYXJhdGUgc2Vjb25kIGZhY3Rvci4gCgpgYGB7cn0KcDE8LWdncGxvdChkYXRhPXN1YnNldChhLCBjYW5kaWRhPT0iWWVzIiksIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvcj1ISVYpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiLAogICAgICAgdGl0bGU9IkNhbmRpZGEiKSsKICBzY2FsZV9zaGFwZV9tYW51YWwoCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApCnAyPC1nZ3Bsb3QoZGF0YT1zdWJzZXQoYSwgY2FuZGlkYSA9PSJObyIpLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3IgPSBISVYpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiLAogICAgICAgdGl0bGUgPSAiTm8gQ2FuZGlkYSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCgKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcF9jb21iaW5lZCA8LSBwMSArIHAyICsgcGxvdF9sYXlvdXQoZ3VpZGVzPSdjb2xsZWN0JykgJiB0aGVtZV9xaygpICYgc2NhbGVfY29sb3JfdWNoaWNhZ28oKQpwX3Jhbmdlc194IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzJdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UpCgpwX3Jhbmdlc195IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc195W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgICAgICAgICAgICAgICBnZ3Bsb3RfYnVpbGQocF9jb21iaW5lZFtbMl1dKSRsYXlvdXQkcGFuZWxfc2NhbGVzX3lbWzFdXSRyYW5nZSRyYW5nZSkKcF9jb21iaW5lZCAmIAogIHhsaW0obWluKHBfcmFuZ2VzX3gpLCBtYXgocF9yYW5nZXNfeCkpICYgCiAgeWxpbShtaW4ocF9yYW5nZXNfeSksIG1heChwX3Jhbmdlc195KSkKYGBgCgpgYGB7cn0KcDE8LWdncGxvdChkYXRhPXN1YnNldChhLCBISVY9PSJZZXMiKSwgYWVzKHg9Tk1EUzEsIHk9Tk1EUzIsIGNvbG9yPWNhbmRpZGEpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiLAogICAgICAgdGl0bGU9IkhJVisiKSsKICBzY2FsZV9zaGFwZV9tYW51YWwoCiAgICAgICAgICAgICAgICAgICAgIGd1aWRlID0KICAgICAgICAgICAgICAgICAgICAgICAgIGd1aWRlX2xlZ2VuZChsYWJlbC50aGVtZSA9IGVsZW1lbnRfdGV4dChzaXplPTYpLCAKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICB0aXRsZT1OVUxMKQogICAgICAgICAgICAgICAgICAgICApCnAyPC1nZ3Bsb3QoZGF0YT1zdWJzZXQoYSwgSElWID09Ik5vIiksIGFlcyh4PU5NRFMxLCB5PU5NRFMyLCBjb2xvciA9IGNhbmRpZGEpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiLAogICAgICAgdGl0bGUgPSAiSElWLSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCgKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcF9jb21iaW5lZCA8LSBwMSArIHAyICsgcGxvdF9sYXlvdXQoZ3VpZGVzPSdjb2xsZWN0JykgJiB0aGVtZV9xaygpICYgc2NhbGVfY29sb3JfdWNoaWNhZ28oKQpwX3Jhbmdlc194IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzJdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UpCgpwX3Jhbmdlc195IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc195W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgICAgICAgICAgICAgICBnZ3Bsb3RfYnVpbGQocF9jb21iaW5lZFtbMl1dKSRsYXlvdXQkcGFuZWxfc2NhbGVzX3lbWzFdXSRyYW5nZSRyYW5nZSkKcF9jb21iaW5lZCAmIAogIHhsaW0obWluKHBfcmFuZ2VzX3gpLCBtYXgocF9yYW5nZXNfeCkpICYgCiAgeWxpbShtaW4ocF9yYW5nZXNfeSksIG1heChwX3Jhbmdlc195KSkKYGBgCgpgYGB7cn0KcDE8LWdncGxvdChkYXRhPXN1YnNldChhLCBISVY9PSJZZXMiKSwgYWVzKHg9Tk1EUzEsIHk9Tk1EUzIsIGNvbG9yID0gc21va2luZ19jdXJyZW50KSApKwogIGdlb21fcG9pbnQoKSsKICBzdGF0X2VsbGlwc2UoKSsKICBsYWJzKHk9Ik1EUzIiLCB4PSJNRFMxIiwKICAgICAgIHRpdGxlPSJISVYrIikrCiAgc2NhbGVfc2hhcGVfbWFudWFsKAogICAgICAgICAgICAgICAgICAgICBndWlkZSA9CiAgICAgICAgICAgICAgICAgICAgICAgICBndWlkZV9sZWdlbmQobGFiZWwudGhlbWUgPSBlbGVtZW50X3RleHQoc2l6ZT02KSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdGl0bGU9TlVMTCkKICAgICAgICAgICAgICAgICAgICAgKQpwMjwtZ2dwbG90KGRhdGE9c3Vic2V0KGEsIEhJViA9PSJObyIpLCBhZXMoeD1OTURTMSwgeT1OTURTMiwgY29sb3IgPSBzbW9raW5nX2N1cnJlbnQpICkrCiAgZ2VvbV9wb2ludCgpKwogIHN0YXRfZWxsaXBzZSgpKwogIGxhYnMoeT0iTURTMiIsIHg9Ik1EUzEiLAogICAgICAgdGl0bGUgPSAiSElWLSIpKwogIHNjYWxlX3NoYXBlX21hbnVhbCgKICAgICAgICAgICAgICAgICAgICAgZ3VpZGUgPQogICAgICAgICAgICAgICAgICAgICAgICAgZ3VpZGVfbGVnZW5kKGxhYmVsLnRoZW1lID0gZWxlbWVudF90ZXh0KHNpemU9NiksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHRpdGxlPU5VTEwpCiAgICAgICAgICAgICAgICAgICAgICkKcF9jb21iaW5lZCA8LSBwMSArIHAyICsgcGxvdF9sYXlvdXQoZ3VpZGVzPSdjb2xsZWN0JykgJiB0aGVtZV9xaygpICYgc2NhbGVfY29sb3JfdWNoaWNhZ28oKQpwX3Jhbmdlc194IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzJdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc194W1sxXV0kcmFuZ2UkcmFuZ2UpCgpwX3Jhbmdlc195IDwtIGMoZ2dwbG90X2J1aWxkKHBfY29tYmluZWRbWzFdXSkkbGF5b3V0JHBhbmVsX3NjYWxlc195W1sxXV0kcmFuZ2UkcmFuZ2UsCiAgICAgICAgICAgICAgICBnZ3Bsb3RfYnVpbGQocF9jb21iaW5lZFtbMl1dKSRsYXlvdXQkcGFuZWxfc2NhbGVzX3lbWzFdXSRyYW5nZSRyYW5nZSkKcF9jb21iaW5lZCAmIAogIHhsaW0obWluKHBfcmFuZ2VzX3gpLCBtYXgocF9yYW5nZXNfeCkpICYgCiAgeWxpbShtaW4ocF9yYW5nZXNfeSksIG1heChwX3Jhbmdlc195KSkKYGBgCgoKIyMgUmVmZXJlbmNlcwo=